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In  section  1.  we  present  i  study  which  addresses  the  velocity  structure  of  the  crust  and  upper  mantle  in  Southern 
China  with  special  emphasis  on  the  Tibet  region.  Velocity  constraints  on  the  crust  are  derived  from  modeling  regional 
Love  waves  and  Pjj  waves  in  the  time  domain  with  observed  waveforms.  An  average  crustal  thickness  of  70  km  is  obtained 
beneath  the  Tibetan  Plateau  with  a  modest  increase  of  velocity  with  depth.  The  lithospheric  and  upper-mantle  structure  is 
deduced  from  modeling  S  and  SS  triplication  waveform  data  and  relative  travel  times.  S-SS  seismograms  chosen  with 
bounce-points  directly  under  Tibet  allow  remote  sensing  of  this  inaccessible  region.  We  conclude  that  the  upper-mantle 
structure  in  the  entire  region  is  basically  shield-like  below  200  km  (SNA).  However,  the  velocity  of  the  lithosphere  is 
abnormally  slow,  roughly  5%  beneath  TibeL  Modeling  two  station  P^  waveform  data  conforms  the  crustal  and  upper 
mantle  model  The  model  for  Tibet  derived  does  not  have  a  distinct  lid,  and  has  a  positive  velocity  gradient  in  the  crust, 
suggesting  crustal  shortening. 

In  section  2,  broad-band  modeling  along  a  regional  shield  path  was  done  for  the  Harvard  recording  of  the  Saguenay 
earthquake.  The  Saguenay  earthquake,  November  25,  1988,  is  one  of  the  first  large  shield  type  events  recorded  by  a 
broadband  -  high  dynamic  range  instrument,  the  Streckeisen  system,  installed  at  Harvard  station  (HRV).  The  event  is 
sufficiently  large  to  be  well  recorded  teleseismically  and  thus  the  source  characteristics  can  be  determined  by  independent 
means  and  considered  known.  This  allows  a  detailed  study  of  the  propagational  effects  along  this  path,  at  an  epicentral 
distance  of  625  km,  where  the  strengths  of  the  surface  waves  can  be  compared  with  the  crustal  body  phases.  Broadband 
modeling  using  standard  analytical  techniques  and  flat  layaed  models  works  amazingly  well  over  the  period  range  of  0.5  to 
20  seconds.  A  detailed  strategy  for  modeling  broadband  regional  phases  is  given  in  terms  of  a  decomposition  of  the 
vertical  and  radial  seismograms  into  three  segments:  (containing  Pn,  pPn.  sPn.  P^P.  P  coupled  PL  waves);  Sn] 

(containing  Sn.  s5n,  Sm$,  etc.);  and  the  fundamental  Rayleigh  waves.  Information  about  the  upper  crust  is  obtained  from 
the  fundamental  Rayleigh  waves  while  cmstal  thickness  and  velocity  gradients  in  the  mantle  are  obtained  from  Pn\  and  5^. 
This  Particular  crustal  model  has  a  thickness  of  35  km  with  a  sharp  moho  and  a  substantial  gradient  in  the  top  20  km  of  the 
mantle,  0.01  km/sec  per  km  for  both  P  and  S  velocities.  The  mantle  velocities.  a=8.2  and  P=4.55  km/sec  are  slower  than 
expected  for  a  shield  environment  Attenuation  is  not  required  for  waveform  modeling  or  for  absolute  amplitude  estimation. 

In  section  3.  three  component  broad-band  waveforms  of  two  small  earthquakes  near  Upland.  California  recorded  on 
the  Pasadena  broad  band  high  dynamic  range  instrument  were  forward  modeled  to  obtain  useful  Green's  functions.  The 
sensitivity  of  the  synthetic  seismograms  to  differing  layered  structure,  boundary  sharpness,  and  two-dimensional  structures 
was  also  investigated.  We  assumed  that  the  sources  of  these  events  were  both  simple  and  known,  as  determined  from  the 
Caltech-USGS  array  first-motions.  A  trapezoidal  time  function  was  chosen  such  that  the  width  of  the  direct  S-wave  was  well 
modeled.  The  waveforms  were  forward  modeled  using  Generalized  Rays.  Reflectivity,  and  Finite-Difference  techniques.  In 
addition,  estimates  of  moment,  fault  dimension,  and  stress  drop  were  computed. 

The  results  of  the  modeling  exercise  indicate  that  a  simple  layer  over  a  half-space  model  is  a  adequate 
approximation  of  the  upper  crust  along  this  profile.  In  particular,  the  waveforms  are  controlled  by  a  relatively  slow 
surficial  layer.  The  sensitivity  analyses  indicate  that  the  details  of  boundary  sharpness,  and  deep  crustal  structure  are  not 
very  important  at  frequencies  less  than  1  Hz,  and  at  ranges  less  than  50  km.  The  thickness  of  the  surface  layer  was 
constrained  to  be  greater  than  3  km  and  less  than  5  km  thick.  The  boundary  sharpness  sensitivity  study  indicated  that  the 
contact  between  the  surficial  material,  and  underlying  materia]  can  be  smoothed  to  2  km  without  adversely  affecting  the 
synthetic  waveform  shapes.  A  number  of  two-dimensional  finite-difference  calculations  were  performed,  and  it  was  found 
that  a  ridge  structure  beneath  the  recorder,  acted  as  a  lowpass  filter,  and  the  lower  frequency  phases  were  largely  unaffected. 
Other  two-dimensional  models  with  ridges  between  the  source  to  receiver  clearly  did  not  fit  the  data  however,  the  extended 
duration  of  these  synthetic  waveforms  may  be  important  in  modeling  more  complicated  waveforms  observed  in  data  from 
events  south  of  this  study.  Synthetic  seismograms  computed  for  the  best  fitting  model  were  used  to  estimate  a  long  period 
moment  of  6x10^2  dyne-cm  (Ml=4-6)  and  1x10^2  dyne -cm  (Ml=3.7)  with  identical  triangular  source  time  durations  of  0.3 
seconds.  Assuming  the  same  fault  dimension  of  0.4  km  from  standard  scaling,  stress  drop  estimates  of  410  and  70  bars  are 
obtained.  Generally,  we  found  that  it  is  possible  to  model  local  waveforms  to  frequencies  of  1  Hz.,  without  a  complete 
understanding  of  fine  structural  detail.  Resulting  Green's  functions  can  be  useful  in  studying  historic  events,  and  in 
simulation  of  large  events  from  a  given  source  region. 

In  section  4,  we  present  analytic  and  numerical  wave  fields  from  an  off-center  explosion  in  an  embedded  solid 
sphere  This  study  investigates  the  effects  of  explosions  in  asymmetric  source  regions  on  the  excitation  of  seismic  body 
waves  We  give  an  analytic  formulation  for  determining  the  wave  fields  from  an  off-center  explosion  in  an  embedded  solid 
sphere  in  an  elastic  whole-space.  As  expected,  this  geometry  generates  shear  as  well  as  compressional  body  waves.  The 
calculated  wave  fields  show  that  the  SV  and  SH  wave  generation  is  determined  by  the  asymmetry  of  the  source  region.  The 
results  are  compared  with  the  known  analytic  solutions  of  an  explosion  in  an  elastic  whole-space  and  at  the  center  of  an 
elastic  sphere  embedded  in  the  whole-space.  The  radiation  patterns  at  different  periods  for  different  parameters  of  the  media 
suggest  that  the  asymmetry  of  the  source  region  has  significant  effects  on  shorter  period,  but  has  only  minor  effects  on 
long  periods  The  long  period  F  to  S  wave  maximum  amplitude  results  are  in  agreement  with  that  for  explosions  in 
axisymmetric  cavities. 
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Summary 


The  research  performed  under  the  contract,  during  the  period  24  May  1989 
through  31  July  1990,  can  be  divided  into  four  main  topics;  modeling  the  shear- 
velocity  structure  of  the  crust  and  upper  mantle  beneath  Tibetan  Plateau  and 
Southeastern  China  using  Love  waves  and  S  and  SS  arrivals,  broadband  modeling 
along  a  regional  shield  path  where  the  strengths  of  the  surface  waves  can  be 
compared  with  the  crustal  body  phases,  broad-band  modeling  of  local  earthquakes, 
and  the  wave  fields  from  an  off-center  explosion  in  a  embedded  sphere. 

In  section  1 ,  we  present  a  study  which  addresses  the  velocity  structure  of 
the  crust  and  upper  mantle  in  Southern  China  with  special  emphasis  on  the  Tibet 
region.  Velocity  constraints  on  the  crust  are  derived  from  modeling  regional  Love 
waves  and  Pni  waves  in  the  time  domain  with  observed  waveforms.  An  average 
crustal  thickness  of  70  km  is  obtained  beneath  the  Tibetan  Plateau  with  a  modest 
increase  of  velocity  with  depth.  The  lithospheric  and  upper-mantle  structure  is 
deduced  from  modeling  S  and  SS  triplication  waveform  data  and  relative  travel 
times.  S-SS  seismograms  chosen  with  bounce-points  directly  under  Tibet  allow 
remote  sensing  of  this  inaccessible  region.  We  conclude  that  the  upper-mantle 
structure  in  the  entire  region  is  basically  shield-like  below  200  km  (SNA). 

However,  the  velocity  of  the  lithosphere  is  abnormally  slow,  roughly  5%  beneath 
Tibet.  Modeling  two  station  Pn\  waveform  data  conforms  the  crustal  and  upper 
mantle  model  The  model  for  Tibet  derived  does  not  have  a  distinct  lid,  and  has  a 
positive  velocity  gradient  in  the  crust,  suggesting  crustal  shortening. 

In  section  2,  broad-band  modeling  along  a  regional  shield  path  was  done  for 
the  Harvard  recording  of  the  Saguenay  earthquake.  The  Saguenay  earthquake, 
November  25,  1988,  is  one  of  the  first  large  shield  type  events  recorded  by  a 
broadband  -  high  dynamic  range  instrument,  the  Streckeisen  system,  installed  at 
Harvard  station  (HRV).  The  event  is  sufficiently  large  to  be  well  recorded 
teleseismically  and  thus  the  source  characteristics  can  be  determined  by  independent 
means  and  considered  known.  This  allows  a  detailed  study  of  the  propagational 
effects  along  this  path,  at  an  epicentral  distance  of  625  km,  where  the  strengths  of 
the  surface  waves  can  be  compared  with  the  crustal  body  phases.  Broadband 
modeling  using  standard  analytical  techniques  and  flat  layered  models  works 
amazingly  well  over  the  period  range  of  0.5  to  20  seconds.  A  detailed  strategy  for 
modeling  broadband  regional  phases  is  given  in  terms  of  a  decomposition  of  the 
vertical  and  radial  seismograms  into  three  segments:  Pnl  (containing  Pn,pPn,  sPn, 
PmP,  P  coupled  PL  waves);  Sni  (containing  Sn,  sSn,  SmS,  etc.);  and  the 
fundamental  Rayleigh  waves.  Information  about  the  upper  crust  is  obtained  from 
the  fundamental  Rayleigh  waves  while  crustal  thickness  and  velocity  gradients  in 
the  mantle  are  obtained  from  Pn\  and  5ni-  This  Particular  crustal  model  has  a 
thickness  of  35  km  with  a  sharp  moho  and  a  substantial  gradient  in  the  top  20  km  of 
the  mantle,  0.01  km/sec  per  km  for  both  P  and  S  velocities.  The  mantle  velocities, 
a=8.2  and  p=4.55  km/sec  are  slower  than  expected  for  a  shield  environment. 
Attenuation  is  not  required  for  waveform  modeling  or  for  absolute  amplitude 
estimation. 

In  section  3,  three  component  broad-band  waveforms  of  two  small 
earthquakes  near  Upland,  California  recorded  on  the  Pasadena  broad-band  high 
dynamic  range  instrument,  were  forward  modeled  to  obtain  useful  Green's 
functions.  The  sensitivity  of  the  synthetic  seismograms  to  differing  layered 
structure,  boundary  sharpness,  and  two-dimensional  structures  was  also 
investigated.  We  assumed  that  the  sources  of  these  events  were  both  simple  and 
known,  as  determined  from  the  Caltech-USGS  array  first-motions.  A  trapezoidal 
time  function  was  chosen  such  that  the  width  of  the  direct  S-wave  was  well 
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modeled.  The  waveforms  were  forward  modeled  using  Generalized  Rays, 
Reflectivity,  and  Finite-Difference  techniques.  In  addition,  estimates  of  moment, 
fault  dimension,  and  stress  drop  were  computed. 

The  results  of  the  modeling  exercise  indicate  that  a  simple  layer  over  a  half¬ 
space  model  is  a  adequate  approximation  of  the  upper  crust  along  this  profile.  In 
particular,  the  waveforms  are  controlled  by  a  relatively  slow  surficial  layer.  The 
sensitivity  analyses  indicate  that  the  details  of  boundary  sharpness,  and  deep  crustal 
structure  are  not  very  important  at  frequencies  less  than  1  Hz,  and  at  ranges  less 
than  50  km.  The  thickness  of  the  surface  layer  was  constrained  to  be  greater  than  3 
km  and  less  than  5  km  thick.  The  boundary  sharpness  sensitivity  study  indicated 
that  the  contact  between  the  surficial  material,  and  underlying  material  can  be 
smoothed  to  2  km  without  adversely  affecting  the  synthetic  waveform  shapes.  A 
number  of  two-dimensional  finite-difference  calculations  were  performed,  and  it 
was  found  that  a  ridge  structure  beneath  the  recorder,  acted  as  a  lowpass  filter,  and 
the  lower  frequency  phases  were  largely  unaffected.  Other  two-dimensional 
models  with  ridges  between  the  source  to  receiver  clearly  did  not  fit  the  data 
however,  the  extended  duration  of  these  synthetic  waveforms  may  be  important  in 
modeling  more  complicated  waveforms  observed  in  data  from  events  south  of  this 
study.  Synthetic  seismograms  computed  for  the  best  fitting  model  were  used  to 
estimate  a  long  period  moment  of  6X1022  dyne-cm  (Ml=4.6)  and  lxlO22  dyne-cm 
(Ml=3.7)  with  identical  triangular  source  time  durations  of  0.3  seconds.  Assuming 
the  same  fault  dimension  of  0.4  km  from  standard  scaling,  stress  drop  estimates  of 
410  and  70  bars  are  obtained.  Generally,  we  found  that  it  is  possible  to  model  local 
waveforms  to  frequencies  of  1  Hz.,  without  a  complete  understanding  of  fine 
structural  detail.  Resulting  Green's  functions  can  be  useful  in  studying  historic 
events,  and  in  simulation  of  large  events  from  a  given  source  region. 

In  section  4,  we  present  analytic  and  numerical  wave  fields  from  an  off- 
center  explosion  in  an  embedded  solid  sphere.  This  study  investigates  the  effects  of 
explosions  in  asymmetric  source  regions  on  the  excitation  of  seismic  body  waves. 
We  give  an  analytic  formulation  for  determining  the  wave  fields  from  an  off-center 
explosion  in  an  embedded  solid  sphere  in  an  elastic  whole-space.  As  expected,  this 
geometry  generates  shear  as  well  as  compressional  body  waves.  The  calculated 
wave  fields  show  that  the  SV  and  SH  wave  generation  is  determined  by  the 
asymmetry  of  the  source  region.  The  results  are  compared  with  the  known  analytic 
solutions  of  an  explosion  in  an  elastic  whole-space  and  at  the  center  of  an  elastic 
sphere  embedded  in  the  whole-space.  The  radiation  patterns  at  different  periods  for 
different  parameters  of  the  media  suggest  that  the  asymmetry  of  the  source  region 
has  significant  effects  on  shorter  period,  but  has  only  minor  effects  on  long 
periods.  The  long  period  P  to  S  wave  maximum  amplitude  results  are  in  agreement 
with  that  for  explosions  in  axisymmetric  cavities. 
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SECTION  1 

Shear-Velocity  Structure  of  the  Crust  and  Upper  Mantle 
Beneath  Tibetan  Plateau  and  Southeastern  China 
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SUMMARY 

This  study  addresses  the  velocity  structure  of  the  crust  and  upper  mantle 
in  Southern  China  with  special  emphasis  on  the  Tibet  region.  Velocity  constraints 
on  the  crust  are  derived  from  modeling  regional  Love  waves  and  P ^  waves  in  the 
time  domain  with  observed  waveforms.  An  average  crustal  thickness  of  70  km  is 
obtained  beneath  the  Tibetan  Plateau  with  a  modest  increase  of  velocity  with 
depth.  The  lithospheric  and  upper-mantle  structure  is  deduced  from  modeling  S 
and  SS  triplication  waveform  data  and  relative  travel  times.  S-SS  sei  .mograms 
chosen  with  bounce-points  directly  under  Tibet  allow  remote  sensing  of  this  inac¬ 
cessible  region.  We  conclude  that  the  upper-mantle  structure  in  the  entire  region  is 
basically  shield-like  below  200  km  (  SNA  ).  However,  the  velocity  of  the  litho¬ 
sphere  is  abnormally  slow,  roughly  5%  beneath  Tibet.  Modeling  two  station  Pn, 
waveform  data  conforms  the  crustal  and  upper  mantle  model.  The  model  for  Ti¬ 
bet  derived  does  not  have  a  distinct  lid,  and  has  a  positive  velocity  gradient  in  the 
crust,  suggesting  crustal  shortening. 

Key  words:  Tibetan  Plateau,  Southeastern  China,  upper  mantle,  shear  velocity 
structure 


INTRODUCTION 

China  is  a  part  of  the  Eurasian  plate,  but  the  margins  of  the  Indian  and 
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Philippine  plates  are  involved  in  the  Himalayas  and  in  the  coastal  ranges  of 
Taiwan,  respectively.  Relative  to  the  north  of  China,  the  Indian  plate  is  moving 
NNE,  and  the  Philippine  plate  is  moving  NW  (  Minster  et  al,  1974  ).  The  arrows 
in  Figure  1  show  the  directions  of  the  motions  of  the  plates  surrounding  China  ( 
Wang  and  Xu,  1985  ).  China  is  composed  mainly  of  four  major  tectonic  provinces, 
the  Tarim  Basin,  Sino-Korea  Craton,  Yangtze  Craton  and  the  Tibetan  Plateau. 
The  Tarim  basin,  Sino-Korea  craton,  and  northern  part  of  Yangtze  craton  are 
precambrian  massives  (  Yang  et  al,  1986;  Zhang  et  al,  1984  ).  The  provinces  are 
characterized  by  large  variations  in  crustal  thickness  as  indicated  in  the  isopach 
map  presented  in  Figure  2  .  The  thickness  of  the  crust  of  the  eastern  part  of  Chi¬ 
na  is  about  35  km,  that  of  northern  China  is  45  km,  and  that  of  the  Tibetan  Pla¬ 
teau  is  more  than  50  km.  Thus,  we  would  expect  the  crustal  and  upper  mantle 
velocity  distributions  beneath  these  four  tectonic  provinces  to  be  much  different. 

The  dynamic  processes  that  lead  to  the  formation  and  maintenance  of  a 
mean  elevation  of  5  km  over  the  Tibetan  Plateau  are  not  well  understood  and  are 
controversial.  Because  of  its  obvious  importance  for  constraining  these  deep 
processes,  the  seismic  velocity  structure  of  the  Tibetan  Plateau  has  been  the  sub¬ 
ject  of  many  studies. 

The  three-year  SinoFranco  joint  research  program  (  1980-1982  )  contribut¬ 
ed  significantly  to  both  the  geology  and  geophysics  of  the  Tibetan  Plateau.  P- 
wave  velocity  profile  of  a  600  km  line  (  from  (28.5t7V,89.0o£')  to  (32. 2° N, 91. 7° E)  ) 
was  obtained  from  seismic  soundings.  The  thickness  of  the  crust  changes  from 
about  40  km  to  75  km  along  this  profile.  The  results  of  the  program  are  pub¬ 
lished  in  a  series  of  papers  (  see  for  example  Him  et  al,  1984;  Teng,  1987  ). 

Some  studies  of  the  seismic  velocity  structure  of  the  Tibetan  Plateau  used 
group  and  phase  velocities  of  fundamental  mode  surface  waves  (  e.  g.  Chen  and 
Molnar,  1981;  Romanowicz,  1982;  1984;  Brandon  and  Romanowicz,  1986;  Feng, 
1982  ).  Other  studies  used  Pn  and  Sn  waves  (  Barazangi  and  Ni,  1982;  Ni  and 
Barazangi,  1983,  Holt  and  Wallace,  1989)  to  derive  the  velocity  of  the  lid  beneath 
the  crust  of  the  plateau.  Recently,  Lyon-Caen  (  1986  )  used  the  travel  times  and 
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waveforms  of  long  period  SH  wave  data  recorded  at  distances  of  10-30  degrees  and 
some  SS-S  waveforms  to  constrain  the  upper  mantle  velocities  down  to  a  depth  of 
400  km  beneath  the  plateau.  She  argues  that  the  Indian  plate  is  not  underthrust¬ 
ing  the  whole  of  the  Tibetan  Plateau  at  the  present  time. 

For  the  other  parts  of  China,  there  are  fewer  studies  reported  in  English. 
For  the  Yangtze  Craton,  there  are  surface  wave  studies  (  Wier,  1982;  Feng,  1982  ). 
Shedlock  and  Roecker  (  1987  )  used  travel  time  inversion  to  study  the  elastic  wave 
velocity  structure  of  the  crust  and  upper  mantle  beneath  the  Sino-Korea  craton. 

Wang  and  Yao  (  1989  )  gave  one  upper  mantle  shear  velocity  model  for  the 
Tibetan  Plateau,  and  one  upper  mantle  shear  velocity  model  for  North  China  ( 
Tarim  basin  and  Sino-Korea  craton  )  by  modeling  the  long-period  SS-S  waveforms 
recorded  at  distances  of  30-60  degrees. 

The  purpose  of  this  study  is  to  constrain  the  upper  mantle  shear  velocity 
structure  of  southern  China  by  using  the  travel  times  and  waveforms  of  Love 
waves  recorded  at  distances  of  7-23  degrees;  long-period  S  waves  recorded  at  dis¬ 
tances  of  15-30  degrees:  and  long-period  SS-S  waves  recorded  at  distances  of  30-60 
degrees.  Recent  work  using  long-period  body  wave  data  have  demonstrated  that 
the  travel  times  and  waveforms  of  SH  waves  recorded  at  distances  up  to  30  de¬ 
grees  and  of  SS-S  waves  recorded  between  30  and  60  degrees  can  provide  some 
constraints  on  the  large-scale  velocity  structure  of  the  upper  670  km  of  the  mantle 
and  particularly  of  the  upper  400  km  (  Grand  and  Helmberger,  1984a;  1984b;  Rial 
et.  al,  1984  ).  In  this  distance  range  the  waveforms  are  controlled  by  the  interfer¬ 
ence  of  phases  whose  turning  points  are  in  the  regions  above  the  400  km  discon¬ 
tinuity,  between  the  400  km  and  670  km  discontinuities,  and  below  the  670  km 
discontinuity.  The  differential  travel  times  of  SS  waves  and  S  waves  are  controlled 
mainly  by  the  velocity  structure  of  the  neighborhood  of  the  bounce  points  of  the 
SS  phases.  This  property  makes  it  possible  to  constrain  the  shear  velocity  struc¬ 
ture  of  an  area  not  big  enough  for  the  pure  path  data.  The  search  of  models  can 
be  achieved  by  comparing  the  synthetics  with  the  data.  The  WKBJ  method  ( 
Chapman,  1978  ),  which  is  discussed  at  length  by  Grand  and  Helmberger  (  1984a 


),  was  used  for  constructing  synthetic  5  and  SS  waves,  and  the  mode  summation 
method  (  Harkrider,  1964;  1970  )  was  used  in  studying  Love  waves. 

DATA  AND  METHODS 

In  this  study,  we  used  travel  times  and  waveforms  of  long-period  SH  wave 
data  from  44  earthquakes,  magnitude  4.5  or  larger,  from  1965  to  1985,  within  or 
around  China  (  Table  1  ),  as  recorded  at  the  HTkSSTV  stations  around  China  ( 
ANP,  BAG,  HKC,  KBL,  LAH,  MAT,  MSH,  NDI,  NIL,  QUE,  SEO,  SHL,  SHK  ) 
and  along  with  some  Chinese  stations  (  GYA,  GZH,  LZH,  XAN  ).  In  Table  1,  the 
source  information  is  from  Bulletin  of  the  International  Seismological  Centre  ( ISC 
),  otherwise  indicated.  The  seismograms  have  high  signal  to  noise  ratios,  and  are 
not  necessarily  from  earthquakes  with  known  source  mechanisms.  There  are  only  a 
limited  number  of  Tibetan  earthquakes  carefully  studied  (e.g.  Ni  and  Barazangi, 
1984;  Molnar  and  Chen,  1983;  Tapponnier  and  Molnar,  1977;  Baranowski  et.  al 
1984).  Because  of  lack  of  local  station  net,  it  is  almost  impossible  to  know  the 
source  mechanisms  of  the  smaller  Tibetan  earthquakes,  which  generated  on-scale 
surface  waves.  However,  we  can  infer  roughly  the  source  mechanism  from  the  S- 
waveforms  (  Langston  and  Helmberger,  1975  ).  This  rough  source  mechanism  is 
sufficient  for  S-,  SS- wave  study  as  demonstrated  by  Grand  and  Helmberger  (  1985 
).  The  source  depth  effects  on  S-  and  SS-waveforms  are  corrected  for  by  using 
teleseismic  S-waveforms  as  the  source  time  history  since  the  ray  parameter 
changes  very  slowly  with  distance  for  distances  greater  than  15.0  degrees.  The 
orientation  of  the  fault  effects  the  amplitude  ratio  of  SS-  to  S-waves,  which  was 
not  used.  These  various  approximations  have  been  used  in  previous  SS-S  studies, 
see  for  example  Rial  et.  al  (  1984  ). 


TABLE  1 


The  crustal  velocity  structure  has  not  been  well  determined  in  the  Tibetan 
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Plateau.  To  constrain  this  portion  of  the  model  we  used  relatively  high-frequency 
Love  waves.  But  instead  of  applying  the  conventional  dispersion  analysis,  we 
chose  to  match  the  Love-wave  synthetics  to  the  observations  directly.  The  ap¬ 
propriateness  of  the  model  is  then  judged  on  the  overall  fit  of  the  synthetics  to  the 
observed  waveform  in  absolute  travel  time. 

Since  Love  wave  data  are  not  normally  used  in  this  fashion,  we  start  with 
a  brief  sensitivity  study  to  test  the  resolving  power  of  this  approach  by  discussing 
some  numerical  experiments.  Figures  3  and  4  summerize  the  results  of  these  nu¬ 
merical  experiments.  Figure  3  shows  the  sensitivity  of  synthetics  to  changes  in  the 
models,  source  depths,  and  source  mechanisms.  Figure  4  displays  the  synthetics 
appropriate  for  different  Moho  discontinuities. 

In  Figure  3,  the  left  column  are  the  models  used  to  generate  the  synthetics 
of  fundamental  Love  waveforms  at  various  depths  and  mechanisms.  These  syn¬ 
thetics  were  generated  at  a  distance  of  1000  km,  assuming  a  the  source  time  histo¬ 
ry  of  a  (1,  1,  1)  trapezoid.  All  synthetics  begin  at  260  seconds.  The  crustal  thick¬ 
ness  is  60  km  for  all  models.  The  average  velocities  of  the  crust  are  same,  3.6 
km/sec.  The  mantle  velocities  are  same,  4.6  km/sec,  except  models  g  which  is  4.4 
km/sec  and  h  which  is  4.8  km/sec.  The  attenuation,  Q0,  used  is  300  for  the  crust, 
and  1000  for  the  mantle.  Comparing  the  synthetics  for  these  various  sensitivities 
we  conclude  that:  the  mantle  velocity  has  very  little  effect  on  both  travel  time  and 
waveforms,  see  the  synthetics  of  modeb  g  and  h;  the  lowest  velocity  of  the  crust 
seems  to  dominate  the  travel  time  of  Love  waves,  see  models  b,  c,  f  and  g;  the 
velocity  gradient  does  not  have  much  effect  on  waveforms;  the  source  depth  has  a 
profound  effect  on  the  waveforms  for  complicated  crustal  modeb,  and  less  effect 
for  simpler  modeb;  source  mechanism  effects  travel  time  only  weakly  compared  to 
shallow  velocity  structure. 

Figure  4  displays  some  properties  of  regional  Love  waves  when  crossing  a 
plateau  boundary  similar  to  the  expected  geometry  of  the  Tibetan  Plateau.  The 
column  on  the  left  displays  the  source-receiver  geometry  in  a  simple  idealized  crust 
where  all  the  modeb  have  the  velocities  given  at  the  bottom.  The  dbtance  is  set 
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at  1000  km,  and  strike-slip  source  is  assumed  at  a  depth  of  9.8  km.  These  syn¬ 
thetics  were  generated  with  generalized  ray  theory  and  finite-difference  method,  see 
Helmberger  and  Vidale  (  1989  ).  The  Love- waveforms  are  quite  similar  which  lead 
us  to  conclude  that  the  crustal  thickness  and  dipping  Moho  do  not  have  much 
effects  for  these  geometries.  Thus,  regional  Love  waveforms  can  be  used  as  a  con¬ 
straint  on  the  velocity  structure  of  the  upper  crustal  layers  for  paths  crossing  the 
plateau  with  the  above  restrictions  on  geometry. 

RESULTS 

The  upper  mantle  shear  velocity  structures  we  derived  are  given  in  Table  2 
and  Figure  5.  TIP  is  the  mode!  for  the  Tibetan  Plateau,  and  ECH  that  for  the 
eastern  part  of  China,  part  of  Yangtze  and  Sino-Korean  Cratons. 

Tibetan  Plateau 

It  is  convenient  to  break  this  discussion  into  two  sections,  namely  crustal 
and  upper-mantle,  where  the  former  is  constrained  by  regional  long-period  Love 
waves  and  the  latter  by  5-  and  SS-waveform  data. 


TABLE  2 


Crust 

Our  preferred  velocity  model  is  given  in  Table  2  with  a  3.75  km  surface 
layer  of  2.55  km/sec  (  same  as  that  of  Chen  and  Molnar  (  1981  )  ),  follow  by  a 
16.25  km  of  3.5  km /sec,  20  km  of  3.7  km/sec,  and  30  km  of  3.8  km/sec. 

The  source-station  paths  used  in  this  Love  wave  investigation  are  given  in 
Figure  6,  along  with  two  P ^  paths  discussed  later.  It  is  difficult  to  measure  the 
initial  arrival  time  of  Love  waves,  so  the  first  large  pulse  was  used  to  denote  rela¬ 
tive  timing  and  the  synthetic  has  been  aligned  accordingly.  Figure  6  indicates  the 
number  of  seconds  the  synthetic  is  faster  than  the  observed  for  that  particular  ray 


fi 


path.  The  comparison  of  the  data  and  synthetics  is  given  in  Figures  7,  and  8. 

In  Figure  7,  we  display  the  data,  whose  sources  are  known  (  Holt  and  Wal¬ 
lace,  1989  ),  along  with  the  synthetics  of  proposed  models,  where  TIP  indicates 
the  synthetics  appropriate  for  the  model  developed  in  this  study,  CHEN  associat¬ 
ed  with  model  S7  of  Chen  and  Molnar  (  1981  ),  and  CHUN  computed  from  the 
model  derived  by  Chun  and  McEvilly  (  1985  ).  The  data  and  synthetics  are 
displayed  in  absolute  travel  time.  From  this  figure,  the  waveform  matches  of  TIP 
are  better  than  the  others,  and  the  travel  times  at  stations  NDI  and  SHL  are 
better.  However,  for  the  paths  of  event  24  to  NIL,  and  event  13  to  LAH,  the  trav¬ 
el  times  of  TIP  are  13  seconds  faster  than  the  data  and  are  not  as  good  as  CHEN. 
CHUN  is  slow.  In  Figure  8,  the  time  shift  was  applied  for  each  comparison.  In 
this  figure,  the  depth  used  to  generate  the  synthetic  is  given.  The  depths  may  be 
different  from  the  corresponding  ones  given  in  Table  1.  Because  the  source  depth 
is  poorly  known,  we  just  compared  the  data  of  unknown  source  depth  with  the 
data  of  known  source  depth  to  estimate  the  source  depth.  The  basic  assumption  is 
that,  the  source  depth  should  be  roughly  the  same  if  the  two  data  are  roughly 
alike.  A  strike-slip  source  with  a  time  history  (  2,  2,  2  )  trapezoid  was  used  to  gen¬ 
erate  all  synthetics.  From  these  two  figures,  the  travel  times  and  the  waveforms  fit 
the  data  reasonably  well,  although  TIP  is  faster,  or  slower  for  certain  paths.  We 
consider  this  model  to  be  an  average  velocity  distribution  of  the  crust  of  the  Ti¬ 
betan  Plateau  which  is  considerably  faster  than  previous  proposed  models. 

Upper-mantle 

The  upper-mantle  shear  velocity  distribution  of  the  Tibetan  Plateau 
beneath  the  crust  is  constrained  by  11  S-waveforms  with  distances  from  15.5°  to 
27.2°,  and  16  SS-S  waveforms  with  distances  from  35.3°  to  59.7°  .  Some  of  the  S- 
wave  data  are  from  Lyon-Caen’s  paper  (  1986  ),  and  some  of  the  SS-S  wave  data 
are  from  Wang  and  Yao’s  paper  (  1989  ).  The  travel  times  of  the  S-waves  are  not 
used,  for  they  are  rather  scattered,  and  a  large  portion  of  their  paths  lie  outside 
the  Tibetan  Plateau.  Note  that  the  SS-S  waveform  data  is  controlled  by  triplica- 
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tion  positions  or  differential  times  wbt-e  rays  share  common  paths  near  the  source 
and  receiver  and  thus  much  less  suspectable  to  travel  time  offsets  caused  by  lateral 
variation  (  Grand  and  Helmberger  1984a  ).  The  ray  paths  of  the  S-wave  data  and 
the  bounce  points  of  the  SS- wave  cata  are  giveu  in  Figure  0a.  The  triangles 
represent  the  stations;  stars  the  events  generating  ihe  S-waves;  circles  the  bounce 
points  of  the  SS-waves  for  distances  greater  than  45  degrees;  and  squares  the 
bounce  points  of  the  SS-waves  for  distances  less  than  45  degrees.  In  order  to  avoid 
the  strong  lateral  heterogeneity,  we  did  not  use  the  SS-wave  data  with  bounce 
points  near  the  boundary  of  the  Tibetan  Plateau. 

The  starting  model  for  the  upper  mantle  shear  velocity  structure  of  the  Ti¬ 
betan  Plateau  is  TNA  (  Tectonic  North  America,  Grand  and  Helmberger,  1984a  ), 
with  the  crustal  model  derived  above.  After  comparisons  of  synthetics  and  the 
data  for  several  dozens  of  models,  following  a  trial-and-error  procedure,  we  ob¬ 
tained  model  TIP,  see  Figure  5  and  Table  2.  This  model  has  a  60  km  thick  lid  of 
4.6  km/sec,  and  40  km  of  4.55  km/sec  below.  It  does  not  have  a  distinct  low  velo¬ 
city  zone,  and  is  shield-like  below. 

Figure  10  shows  the  comparison  of  the  observed  waveforms  with  synthet¬ 
ics.  Because  most  of  the  5-wave  data  do  not  correspond  to  pure  paths  and  sample 
an  area  with  large  heterogeneities,  it  is  difficult  to  say  which  part  of  the  data 
gives  the  most  definitive  information  about  the  Tibetan  Plateau.  However,  by 
conpparing  the  paths  indicated  in  Figure  0a  we  see  that  the  shallowest  structure  is 
sampled  best  by  the  western  station,  NIL.  Almost  all  of  the  data  and  synthetics 
have  three  arrivals  that  are  due  to  three  sets  of  rays,  one  bottoming  above  the  405 
km  discontinuity  (  branch  AB  ),  one  bottoming  below  the  405  km  discontinuity 
and  above  the  670  km  discontinuity  (  branches  BC  and  CD  ),  and  another  bot¬ 
toming  below  the  670  km  discontinuity  (  branches  DE  and  EF).  We  will  generally 
refer  to  branches  by  single  letters,  which  will  denote  the  two  arrivals  forming  the 
cusp  designated  by  the  letter  (  Figure  11  ).  With  long-period  data  the  two 
geometric  arrivals  forming  a  branch  usually  cannot  be  distinguished.  When  there 
is  possible  confusion  in  the  branch  name,  we  will  indicate  a  full  designation.  A 
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good  fit  of  the  waveforms  and,  therefore,  differential  travel  times  of  branches  A, 
B,  C,  D,  E  and  F  to  the  data  at  distances  of  17.2°  ,  20.1°  ,  24.0°  ,  24.6®  ,  25.3®,  and 
26.7®  is  achieved.  Note  that  the  A  branch  is  missing  from  the  LZH  data  (  18.4°  ) 
indicating  a  shallow  heterogeneity  along  the  path  (  Grand  and  Helmberger,  1985  ). 
The  differential  travel  time  between  the  branches  A  and  C  of  the  data  recorded  at 
NIL  at  a  distance  of  15.5°  is  3.5  seconds  faster  than  that  of  the  synthetics,  and 
that  of  the  data  also  recorded  at  NIL  at  a  distance  of  16.4°  is  2  seconds  faster  than 
that  of  synthetics.  Perhaps  these  differences  were  caused  by  the  local  faster  struc¬ 
ture  and  thinner  crust  thickness  near  the  station  NIL.  The  amplitude  of  the  E 
branch  for  QUE  data  at  a  distance  of  21.0°  could  be  due  to  the  missing  A  branch. 

The  travel  times  and  waveforms  of  S- waves  alone  do  not  yield  definitive 
constraints  on  the  velocity  distribution  of  the  Tibetan  Plateau,  because  a  large 
part  of  the  ray  paths  are  outside  the  region.  Figure  12  shows  the  SS-S  vrave  data 
with  synthetics  of  the  model  TIP.  The  fits  of  the  waveforms  and  differential  trav¬ 
el  times  between  the  first  and  the  second  arrivals  (  F  and  D  branch  ),  and  the  S- 
wave  of  the  data  recorded  at  SEO  at  a  distance  of  54.5  degrees  with  the  synthetics 
indicate  that  the  shear  velocity  distribution  beneath  the  405  km  discontinuity  is 
the  same  as  that  beneath  the  Canadian  Shield  (  Grand  and  Helmberger,  1984  ). 
The  source  functions  used  for  the  synthetics  of  SS  waveforms  are  the  teleseismic  5 
waveforms.  The  advantage  of  using  the  teleseismic  S-waveform  as  a  source  func¬ 
tion  is  that  the  importance  of  the  source  mechanisms  and  source  depth  is  almost 
completely  eliminated.  The  waveforms  of  5-waves  are  the  sum  of  direct  5-  and 
«S-waves.  The  travel  time  difference  between  these  two  phases  is  controlled  by  the 
source  depth  and  the  ray  parameter,  or  take-off  angle.  The  amplitude  ratio  is 
determined  by  the  source  mechanism  (  Langston  and  Helmberger,  1975  ).  The 
different  S-waveforms  at  different  distances  are  due  to  the  ray  parameters.  The 
ray  parameter  changes  very  slowly  with  distance  beyond  15  degrees,  becoming 
nearly  a  constant  for  distances  beyond  25  degrees.  Thus  the  effects  of  source  depth 
and  source  mechanism  are  almost  the  same  on  triplication  S-waves  as  on  telese¬ 
ismic  S-  waves. 
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The  differential  travel  time  of  the  C  branch  and  S-wave  of  the  data  record¬ 
ed  at  MSH  at  a  distance  of  35.3°  is  about  5  seconds  faster  than  that  of  the  syn¬ 
thetics,  and  the  branch  A  is  absent.  The  bounce  point  for  it  is  the  far  west  point 
indicated  in  Figure  9a.  The  ray  of  the  C  branch  in  this  case  travels  almost  entire¬ 
ly  in  the  upper  405  km  with  only  a  very  small  part  of  it  going  below  the  405  km 
discontinuity.  A  large  part  of  the  time  difference  is  due  to  the  more  than  10  de¬ 
grees  of  horizontal  travel  distance  in  a  shield-like  region  (  Rial  et.  al,  1984  )  out¬ 
side  of  the  Tibetan  Plateau.  Also,  the  lateral  heterogeneity  could  cause  some 
difference  in  the  differential  travel  time  of  SS-  and  S-waves.  The  differential  travel 
times  and  waveforms  of  SS-  and  S-waves  of  the  HKC  (  35.7°  )  data,  see  Figure  12, 
are  matched  very  nicely  for  the  branches  C  and  E,  although  the  A  branch  is  not 
clear.  The  beginning  waveform  of  the  HKC  station  at  a  distance  of  52.0°  is  most 
likely  S.  If  it  is  a  combination  of  S  and  sS,  the  source  depth  should  be  about  80 
km,  and  the  amplitude  of  the  second  one  would  not  likely  be  so  big.  If  we  use  all 
of  the  beginning  waveform  as  the  S-wave,  or  source  function,  the  synthetics  for 
the  SS-wave  are  very  good,  but  we  could  not  label  the  branches.  The  waveforms 
of  the  SS- waves  of  HKC  at  a  distance  of  51.0°  do  not  agree  well  with  that  of  the 
synthetics,  but  the  match  of  the  differential  travel  time  is  not  bad.  The  arrival 
between  the  S-  and  SS- waves  on  the  data  HKC  51.0°,  GZH  51.8°,  HKC  52  0°,  SEO 
54.5°,  and  HKC  59.7°  is  the  phase  ScS. 

The  synthetic  fits  to  this  dataset  proved  particularly  difficult.  The  reason  is 
in  the  inherent  complexity  of  the  structure  itself  in  terms  of  heterogeneity. 
Nevertheless,  it  appears  that  the  derived  model  for  the  structure  beneath  the  Ti¬ 
betan  Plateau  fits  the  S-SS  waveform  data  better  than  existing  models,  as  will  be 
discussed  in  more  detail  later. 

Eastern  China 

The  data  coverage  of  the  eastern  portion  of  China  is  presently  lacking  but 
will  become  more  complete  when  the  new  digital  array  becomes  operational.  The 
ray  paths  of  the  data  used  for  the  study  of  Eastern  China  are  displayed  in  Figure 


12 


13.  It  is  clear  from  the  figure  that  the  region  sampled  includes  the  Yangtze  Craton 
and  the  southern  half  of  the  Sino-Korea  Craton.  Although  the  region  is  relatively 
complex,  we  will  assume  a  uniform  upper-mantle  structure  as  a  first  order  approx¬ 
imation  of  the  area.  The  data  set  we  used  to  constrain  the  upper  mantle  shear 
velocity  distribution  includes  5  S-wave  data  ,  distances  from  19.0°  to  24.5°,  and  10 
55wave  data,  distances  from  33.2°  to  39.7°.  The  model,  ECH,  we  obtained,  is 
given  in  Table  2  and  Figure  5.  The  comparisons  of  the  synthetics  with  the  data 
are  given  in  Figure  14  and  Figure  15. 

In  Figure  14,  we  display  synthetic  and  5  wave  observation  comparisons 
along  with  the  theoretical  responses  and  source  function  estimates.  The  branches 
are  identified  on  the  theoretical  responses  similar  to  Figure  10  discussed  earlier. 
Although  all  the  branches  cross  each  other  in  this  distance  range,  we  can  see 
pulses  of  branches  on  both  the  data  and  the  synthetics.  The  travel  time  of  the 
synthetics  is  3  seconds  faster  than  that  observed  at  HKC  (  19.0°  ),  3  seconds 
slower  than  that  of  ANP  (  20.8°  ),  4  seconds  slower  than  that  of  SHL  (  22.9°  ),  7 
seconds  slower  than  that  of  SEO  (  23.2°  ),  and  9  seconds  slower  than  that  of  ANP 
(  24.5°  )  (  Figure  13  ).  The  fit  of  the  travel  times  is  certainly  reasonable  for  such  a 
big  laterally  heterogeneous  area.  The  synthetics  are  produced  by  convolving  the 
theoretical  responses  (  Earth  responses  )  with  the  effective  source  functions  which 
contain  sS  .  The  source  function  of  ANP  (  20.8°  )  is  the  S-waveform  observed  at 
MAT  (  36.0°,  Figure  15  ).  The  other  source  functions  are  theoretical  predictions 
assuming  a  strike-slip  at  SHL  and  SEO,  and  a  dip-slip  at  ANP  and  HKC.  The 
bottom  two  observations  at  SEO  and  ANP,  showing  the  interference  of  F  and  D 
branches,  are  very  similar  to  those  displayed  in  Figure  17  of  Grand  and  Helm- 
berger  (  1984a  ). 

Figure  15  gives  the  comparisons  of  the  synthetics  of  55waveforms  and  the 
data.  The  55-waveform  recorded  at  a  distance  less  than  45°  gives  not  only  infor¬ 
mation  about  the  bounce  point,  but  also  information  about  a  large  part  of  the 
upper  mantle  structure  along  the  ray  path.  All  of  the  55 wave  data  we  have  used 
are  less  than  40°,  so  this  may  cause  some  difficulties  for  the  large  lateral  hetero- 
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geneity.  Also  it  is  difficult  to  locate  the  E,  F  branches  of  the  model.  The  complexi¬ 
ty  of  the  waveforms  of  SHK  (  35.0°  )  is  due  to  the  depth  of  the  source,  134  km. 
We  have  not  matched  sS-  and  sSS- waves,  so  the  synthetic  is  incomplete.  Note 
that  7  of  the  10  seismograms  used  are  recorded  at  MAT,  and  one  third  of  their 
ray  paths  lie  out  of  continental  China  (  Figure  13  ).  Actually,  for  some  of  the  SS- 
wave  data,  only  one  third  of  their  paths  and  bounce  points  are  in  the  eastern  part 
of  China,  Yangtze  Craton,  and  the  south  half  of  the  Sino-Korea  Craton.  It  can  be 
inferred  from  the  data  and  the  synthetics  (  Figure  14,  and  Figure  15  )  and  the 
travel  time  residuals  (  Figure  13  )  that  the  velocity  of  the  Sino-Korea  Craton  is 
faster  than  that  of  the  Yangtze  Craton,  and  that  the  northern  part  of  Yangtze 
Craton  is  faster  than  the  southern  part.  The  model  we  derived  here  is,  obviously, 
only  a  very  approximate  model  for  the  upper  mantle  shear  velocity  distribution 
for  the  eastern  part  of  China.  Nevertheless,  these  preliminary  results  suggests  that 
the  mantle  beneath  eastern  China  is  predominantly  shield  at  depths  greater  than 
a  few  hundred  kilometers.  The  large  travel  time  off-sets  suggest  strong  shallow  la¬ 
teral  variation.  More  data  will  be  required  to  resolve  the  connection  with  the  sur¬ 
face  geology  and  interplate  interaction. 

DISCUSSION  AND  CONCLUSION 

In  this  section  we  will  briefly  compare  synthetics  associated  with  existing 
models  (  Figure  16  )  proposed  for  the  Tibetan  Plateau  with  some  of  the  key  SS-S 
observations.  The  SS-S  data  is  the  least  contaminated  by  lateral  variation  as  dis¬ 
cussed  earlier  and  thus  the  most  definitive.  This  will  be  followed  by  an  indepth 
discussion  of  the  upper  200  kms  of  the  models  since  this  is  the  region  where 
models  differ  the  most  and  are  the  most  significant  in  terms  of  tectonic  implica¬ 
tions.  The  waveforms  of  SS-S  waves  are  the  results  of  the  interference  of  the  five 
branches  AB,  BC,  CD,  DE  and  EF  (  Figure  11).  Since  the  travel  time  branches  are 
nearly  straight  lines,  we  can  get  a  shear  velocity  model  for  a  homogeneous  struc¬ 
ture  if  we  have  two  ideal  waveforms  of  SS-S  waves  at  a  distance  near  30  and  60 
degrees.  We  present  three  key  data  taken  from  Figure  12,  located  at  roughly  10* 


intervals,  namely,  HKC  at  35.7°,  ANP  at  44.8°,  and  SEO  at  54.5°  in  Figure  17.  The 
waveforms  of  SEO  at  54.5°  are  particularly  meaningful,  because  waveforms  of  the 
different  branchs  have  separated.  At  ANP  the  branchs  are  all  together.  In  order  to 
see  clearly  the  differences  between  the  data  and  synthetics,  we  have  inserted  the 
vertical  lines.  The  bottom  four  traces  of  each  group  are  synthetics  for  proposed 
models  for  this  region. 

Comparing  the  various  models  presented  in  Figure  16  we  note  that  TIP  has 
a  faster  crust  and  a  slower  lid  than  the  other  models.  Thus,  the  Sn  and  Pn  arrival 
time  prediction  for  TIP  are  distinctly  slower  than  others.  Comparing  with  the 
data  at  HKC  (  35.7°  ),  the  branches  C  and  E  of  the  TIP  synthetics  are  about  two 
seconds  faster  than  the  data,  this  means  that  either  the  shear  velocity  of  the 
upper  part,  especially  the  crust,  is  faster  than  the  observed  or  the  crustal  thick¬ 
ness  at  the  bounce  point  is  greater  than  the  model  TIP,  as  discussed  earlier.  The 
C  branch,  and  the  E  branch  of  the  synthetics  of  Lyon-Caen’s  models  (  Ltib  )  are 
7.5,  and  4.5  seconds  slower  than  that  of  the  HKC  data.  The  synthetics  for 
Wang’s  model  are  also  slower  than  the  HKC  data  by  a  few  seconds.  The  branches 
arrive  at  almost  same  time  for  the  ANP  (  44.8°  )  data.  The  synthetics  of  Wang 
and  Ltib  are  slower.  The  F,  D  and  B  branches  of  Ltib  are  3.5,  5  and  16  seconds 
slower  than  those  of  SEO  data  respectively.  The  B  branch  of  Wang  is  11  seconds 
slower  than  that  of  SEO  data.  Lyon-Caen’s  Indian  model  presents  a  comparison 
of  Indian  Shield  with  the  Tibetan  Plateau.  The  slow  B  branch  of  Ltib  and  Wang 
suggest  that  the  velocities  of  the  upper  200  km  of  the  models  are  too  slow. 

To  check  our  results  further  requires  some  absolute  travel  time  constraints 
or  better  knowledge  about  the  earthquake  sources  used  in  terms  of  location  and 
mechanism.  For  instance,  we  used  the  ISC  location  and  origin  time  for  event  22  in 
timing  NDI,  see  Figure  6,  and  found  a  residual  of  +10  seconds.  This  means  that 
the  model  TIP  is  1.5  %  faster  than  the  average  velocity  along  this  path,  but  if  we 
used  the  USGS  location  and  origin  time  for  this  event,  the  travel  time  residual 
would  be  +2  seconds,  which  means  this  model  is  very  good.  We  have  no  particular 
reason  to  say  the  ISC  location  and  origin  time  are  better.  We  choose  to  use  them 
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simply  for  book-keeping.  Since  these  locations  given  by  the  agencies  do  not  use 
depth  phases  in  their  analysis  we  decided  to  do  an  indepth  study  of  one  of  these 
events  to  assess  uncertainties  and  establish  a  few  absolute  travel  time  constraints. 
Event  36  was  selected  for  this  purpose  since  it  also  has  a  few  P ^  waves  on  scale, 
the  origin  time  and  location  used  are  3:2:47.2,  30.506*  N,  88.583*  E,  and  an  epicen- 
tral  depth  of  33  km  assigned  by  NEIS.  This  event  is  big  enough  to  be  recorded 
worldwide,  m»*=5.7,M#=6.2. 

First,  we  derive  the  source  depth  by  modeling  the  teleseismic  P  waveforms. 
The  results  are  given  in  Figure  18  where  the  best  fitting  source  depth  is  10  km. 
We  estimate  the  errors  in  depth  to  be  less  than  3  km,  and  fault  plane  less  than  5°. 
Figure  18  displays  the  best  overall  fit  of  different  runs,  indicating  a  strike  of  160°, 
a  dip  of  60°,  and  a  rake  of  250°  .  Stations  POO  and  CHG  are  within  upper  man¬ 
tle  distances  and  the  P  waveforms  recorded  by  them  are  complicated  by  triplica¬ 
tions,  and  matched  well  by  the  synthetics.  The  model  used  in  these  calculation 
was  derived  by  comparing  TIP  with  a  P-wave  model  for  the  Canadian  Shield, 
namely  S25  (  LeFevre  and  Helmberger,  1989  )  and  SNA  (  Grand  and  Helmberger, 
1984a  ),  and  using  the  relationship  with  shear  velocities  of  rocks  of  basaltic  com¬ 
position  (  Ludwig  et  al,  1970  )  (  Table  3  ).  Some  bad-looking  seismograms  are 
due  to  digitizing  since  the  line-thickness  of  the  recordings  is  3-5  seconds  thick,  and 
the  maximum  amplitudes  of  some  recordings  do  not  exceed  half  of  the  line- 
thickness,  such  as  MUN. 


TABLE  3 


The  results  of  the  relocation,  based  on  the  new  source  depth,  are  given  in 
Table  4.  In  this  table,  STAT  are  station  names;  TPK  are  the  times  we  picked 
from  vertical  component  of  short-period  WWSSN  data  after  minute  mark  in 
second;  TPKP  are  the  times  picked  by  station  operator;  TOBSl  and  TOBS2  are 
the  observed  P  wave  travel  times  after  the  old  and  new  origin  time;  TCOMl  and 
TCOM2  are  the  computed  P  wave  travel  times  for  the  old  epicenter  and  origin 
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time,  and  for  the  new  epicenter  and  origin  time.  DTI  and  DT2  are  the  residual 
travel  times  after  Dziewonski  and  Anderson  (  1983  )  station  correction,  minus  sign 
means  TCOM  is  faster  than  TOBS.  The  standard  deviation  is  1.5  seconds.  This 
error  mainly  comes  from  the  contributions  of  the  three  stations  UME,  HLW  and 
STU,  whose  deviations  are  greater  than  4  seconds.  After  relocation:  the  location 
becomes  30.658®  N,  88.649®  E,  and  the  origin  time  becomes  03:02:43.7.  The  new 
location  is  18  km  from  the  location  given  by  NEIS.  The  origin  time  is  3.5  seconds 
earlier  than  that  given  by  NEIS,  1.3  seconds  earlier  than  that  given  by  ISC.  We 
used  the  modified  P-wave  TIP  model  as  discussed  before  for  the  source  region  and 
JB  for  the  receivers  to  calculate  the  TCOM.  For  a  distance  greater  than  30®,  the 
half  way  travel  time  difference  between  the  model  TIP  and  the  model  JB  is  only 
0.5  seconds.  If  we  used  JB  model  for  the  source  region,  the  new  location  would  be 
30.625®  N  ,  88.632°  E,  only  4  km  from  the  location  given  by  different  source  veloci¬ 
ty  models,  the  origin  time  would  be  03:02:43.7,  only  0.04  seconds  difference.  This 
is  because  the  upper  400  km  of  JB  model  is  slow,  or  the  average  velocity  of  upper 
400  km  of  JB  is  about  same  as  that  of  the  model  TIP. 


As  a  check  on  the  accuracy  of  our  crustal  model,  we  make  synthetics  of  P ^ 
waveforms  and  compare  with  the  data  recorded  at  the  stations  NDI  and  NIL  due 
to  the  Tibetan  earthquake  we  relocated  above.  Figure  19  displays  the  P ^  data 
along  with  the  synthetics  of  the  model  TIP.  The  method  used  tc  make  these  syn¬ 
thetics  is  reflectivity.  The  synthetics  of  NDI  is  shifted  left  1.5  seconds,  that  of 
NIL  is  shifted  left  0.5  seconds.  We  present  only  three  seismograms  here  since  the 
amplitudes  of  the  north  components  of  both  stations  are  very  small  and  the  verti¬ 
cal  component  of  NIL  station  was  not  available.  We  can  see  from  Figure  19  that 
the  fits  are  good.  This  means  that  the  average  velocity  structure  and  the  crustal 
thickness  of  the  Tibetan  Plateau  can  be  approximated  by  the  model  TIP  for  the 
paths  (  Figure  6  ),  although  these  paths  only  sample  the  southern  part  of  the  Pla- 
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teau.  The  time  shifts  given  above  are  largely  caused  by  the  dipping  Moho.  Thin¬ 
ning  the  crust  at  the  receiver  by  20  km  reduces  the  Pn  travel  time  by  2.2  seconds, 
see  Figure  6.  Thus,  considering  the  possible  crustal  thickness  beneath  the  stations, 
and  assuming  that  the  crustal  thickness  beneath  NDI  and  NIL  is  50  km,  we  con¬ 
clude  that  the  model  TIP  is  about  0.7  seconds  faster  for  the  path  to  NDI,  1.7 
seconds  faster  for  that  to  NIL.  This  implies  that  the  P  wave  velocity  of  the  man¬ 
tle  for  these  paths  is  8.23  km/sec  instead  of  the  8.29  km/sec  that  we  used  in  the 
above  flat  layered  Pn,  calculation. 

Barazangi  and  Ni  (  1982  ),  and  Ni  and  Barazangi  (  1983  )  used  Pn-  and  Sn- 
waves  crossing  the  Tibetan  Plateau,  and  obtained  velocities  of  8.42  km/sec  for 
Pn,  and  4.73  km/sec  for  Sn.  They  concluded  that  these  velocities  were  very  simi¬ 
lar  to  those  beneath  the  Himalaya  Mountains  and  the  Indian  shield,  and  suggest¬ 
ed  that  the  Indian  continental  lithosphere  underthrusts  the  Tibetan  Plateau  at  a 
shallow  angle.  However,  our  earlier  experiments  suggest  a  high  velocity  bias  for 
most  recording  geometries  if  one  simply  assumes  a  flat-layered  model.  For  in¬ 
stance,  in  Figure  4,  the  S„  arrives  about  3.5  seconds  earlier  in  the  model  d  than  a 
flat  structure  with  a  crustal  thickness  70  km,  model  f.  A  mantle  velocity  would 
be  4.69  km/sec,  instead  of  4.6  km/sec  used  in  Figure  4.  Thus,  an  average  S„  velo¬ 
city  of  0.1  km/sec  greater  than  the  average  velocity  beneath  the  Tibetan  Plateau 
is  very'  likely  obtained  if  the  only  Sn  data  are  used.  Seismograms  for  two- 
dimensional  models  crossing  this  interesting  region  will  be  discussed  in  a  later  pa¬ 
per. 

Our  results  (  TIP  )  are  consistent  with  that  of  the  attenuation  study  of 
pure  path  long-period  Rayleigh  waves  across  the  Tibetan  Plateau  by  Romanowicz 
(  1984  ).  The  conclusion  of  "no  lid"  for  the  Tibetan  Plateau  is  consistent  with  that 
for  the  central  Chang-Thang  of  Tibet  from  pure  path  phase  velocity  measurement 
of  long  period  Rayleigh  waves  by  Brandon  and  Romanowicz  (  1986  ).  These 
results  are  similar  to  those  obtained  from  modeling  multi- bounce  S-waves  by 
Grand  and  Helmberger  (  1985  ),  who  constructed  a  2D  cross-section  from  Tibet  to 
Europe.  Resolution  from  their  study  begins  at  the  edge  of  Tibet  and  indicates 
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that  the  structure  beneath  the  Tarim  basin  is  similar  to  our  results  for  Tibet, 
namely,  a  shield-like  model  with  a  slow  upper  200  km.  We  do  not  presently  have 
a  detailed  shear  velocity  models  available  for  the  eastern  part  of  China,  except 
that  it  appears  similar  to  other  shields.  Thus,  from  this  study,  and  the  studies  by 
Grand  and  Helmberger  (  1085  ),  and  Rial,  Grand  and  Helmberger  (  1984  )  we  sug¬ 
gest  that  the  upper  mantle  of  the  Eurasian  plate  is  shield-like  below  200  km. 

Although  these  results  are  preliminary  in  nature  it  would  appear  that  the 
deformations  of  the  European  plate  is  confined  mainly  to  the  upper  200  km.  This 
type  of  velocity  distribution  thus  becomes  an  important  piece  of  evidence  to  be 
used  in  deducing  the  tectonics  of  the  area. 

It  is  generally  assumed  that  major  undertrusting  is  occurring  along  the 
Himalayan  arc  as  deduced  from  source  mechanism  studies,  see  Ni  and  Barazangi  ( 
1984  )  and  others.  The  dynamic  processes  that  lead  to  the  formation  and  mainte¬ 
nance  of  a  mean  elevation  of  5  kms  over  some  2  million  square  kms  behind  the 
Himalayan  arc  appears  more  controversial.  Essentially  two  hypothesis  have  been 
put  forward,  namely  underthrusting  and  crustal  shortening.  Underthrusting  in  the 
Pacific  Northwest  results  in  complex  crustal  structures  with  a  large  low  velocity 
zone  associated  with  the  subducted  upper  crustal  section,  see  Langston  (  1977  ). 
We  did  not  find  much  evidence  for  such  structures  in  this  study.  Crustal  shorten¬ 
ing  has  been  suggested  by  Molnar  and  his  associates.  They  argue  that  the  tem¬ 
perature  at  the  base  of  the  crust  is  250  to  300°  higher  than  beneath  platforms,  see 
Chen  and  Molnar  (  1981  ).  Increased  radiogenic  material  per  horizontal  area  is 
thought  to  contribute  to  this  heating.  If  we  compare  the  mantle  shear  velocity  at 
the  top  of  our  model  with  those  of  pure  shield  models  we  obtain  a  reduction  of 
0.1  -  0.2  km/sec,  essentially  4.7  -  4.8  to  4.6  km/sec.  The  temperature  difference  for 
this  velocity  difference  is  300-600°,  assuming  that  the  upper  mantle  is  composed  of 
olivine  (  40  %  ),  clinopyroxene  and  garnet  as  suggested  by  Duffy  and  Anderson  ( 
1989  ).  We  have  (d  0/dT),-*— 3.3X10~Akm /(acc  °C)  for  upper  mantle  materials  at 
1000°Ar  (  Duffy,  personal  communication;  Suzuki  et.  al,  1983  ).  This  temperature 
estimation  are  somewhat  higher  than  the  above  estimates,  see  Molnar  (  1989  ), 


reflecting  the  contrast  in  velocity  models,  see  Figure  16,  where  the  Lyon-Caen’s 
model  has  a  higher  lid  velocity  than  TIP.  Our  results  suggests  that  the  crustal 
structure  is  playing  a  particularly  important  role  in  controlling  the  tectonics  of 
the  Tibetan  Plateau  and  in  surrounding  regions.  Thus  the  higher  resolution  of 
this  complex  structure  is  essential  and  will  be  addressed  in  future  efforts  as  the 
new  digital  observations  from  the  Soviet  Union  and  China  become  available. 
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Table  1.  Earthquakes  and  stations  used  in  this  study. 


Date 

Origin 

Time 

1 

06 

14 

65 

13:17:02.4 

2 

06 

16 

65 

23:49:08.2 

3 

06 

17 

65 

20.14.50.1 

4 

06 

18 

65 

01:18:39.0 

5 

02 

05 

66 

15:12:32.9 

6 

02 

07 

66 

23:06:37.4 

7 

03 

14 

66 

04:42:50.7 

8 

03 

17 

66 

05:44:47.9 

9 

03 

29 

66 

06:12:00.8 

10 

06 

28 

66 

14:00:21.0 

11 

03 

14 

67 

06:58:04.4 

12 

08  15  67 

09:21:03.3 

13 

08  30  67 

11:08:50.0 

14 

09  15  67 

10:32:44.2 

15 

12  22  68 

09:06.35.0 

16 

10  17  69 

01:25:11.5 

17 

06  05  70 

04:53:07.4 

18 

07  30  70 

00:52:20.3 

19 

04  03  71 

04:49  03.1 

20 

05  31  71 

05:13:58.6 

21 

07  22  72 

16:41:02.1 

22 

08  30  72 

15:14:7.5 

23 

08  30  72 

18:47:40.3 

24 

02  07  73 

16.06:25.8 

25 

03  18  75 

18:44:16.3 

26 

05  05  75 

05:18:46.3 

27 

10  03  75 

17:31:35.6 

28 

04  08  76 

02:40:23.9 

29 

05  31  76 

05:08:30.0 

30 

05  31  76 

18:35:05.0 

31 

07  03  78 

08:46:32.0 

32 

03  15  79 

12:52:26.0 

33 

05  20  79 

22:59:11.8 

34 

01  12  80 

15:31:40.0 

35 

02  13  80 

22:09:30.8 

36 

02  22  80 

03:02:43.7 

37 

07  29  80 

12:23:07.7 

38 

01  23  82 

17:37:29.2 

39 

06  15  82 

23.24:28.8 

40 

04  18  83 

10:58:49.1 

41 

12  16  83 

13:15:57.3 

42 

02  01  84 

14:22:09.2 

43 

10  26  84 

20:22:21.8 

44 

10  29  85 

13:13:40.0 

Location 

Depth 

Mag 

(*N  *E  ) 

(  km  ) 

Mb  Ms 

32.09  87.72 

36 

5.2 

32.04  87  .46 

69 

4.8 

32.12  87.76 

15 

5.1 

32.01  87.59 

43 

5.1 

26.22  103.21 

32 

5.6 

30.25  69.89 

28 

5.6 

32.47  97.46 

33 

4.8 

31.60  82.76 

11 

4.7 

37.52  114.99 

S3 

5.3 

27.53  100.08 

12 

5.7 

28.41  94.29 

12 

5.8 

31.05  93.56 

10 

5.5 

31.57  100.31 

8 

4.9 

27.42  91.86 

19 

5.8 

36.25  101.83 

21 

5.5 

23.09  94.70 

124 

6.1 

42.48  78.71 

24 

5.9  6.8 

37.85  55.94 

22 

5.7 

32.16  94.99 

27 

5.6  6.4 

25.22  96.51 

22 

5.2  6.1 

31.38  91.41 

17 

5.4  5.8 

36.65  96.35 

17 

5.5  5.5 

36.56  96.35 

16 

5.5  5.3 

31.50  100.33 

8 

5.9 

35.12  86.54 

31 

5.1  5.8 

33.13  92.84 

8 

5.6  6.1 

30.44  66  41 

24 

5.5  6.4 

40.31  83.72 

10 

6.2  6.4 

24.37  98.62 

25 

5.5  6.2 

24.29  98  68 

20 

5.1  5.5 

32.75  94  08 

8 

4.7  4.9 

23.18  101.09 

6 

5.5  6.2 

29.93  80.27 

16 

5.7  5.7 

33.58  57.26 

14 

5.3  5.8 

36.47  76.86 

74 

6.0  5.4 

30.66  88.65 

10 

6.C  5.4 

29.34  81.21 

3 

5.7  5.5 

31.68  82.28 

25 

6.0  6.5 

31.85  99.92 

7 

27.78  62.07 

44 

6.4  6.3 

39.34  72.96 

35 

5.7  5.9 

34.57  70.48 

44 

5.9  5.9 

39.15  71.35 

7 

5.9  6.2 

36.75  54.81 

13 

6.0  8.0 

Study 

Station 

A 

Phase 

area 

TP 

SHL 

7.5 

Love 

LAH 

11.4 

Love 

TP 

NDI 

9.5 

Love 

TP 

LAH 

11.4 

Love 

TP 

NDI 

9.6 

Love 

EC 

SEO 

23.1 

S 

TP 

ANP 

45.7 

SS 

TP 

SHL 

8.4 

Love 

NDI 

17.8 

Love 

TP 

SHL 

10.0 

Love 

EC 

SHL 

23.0 

S 

TP 

MSH 

35.3 

SS 

EC 

HKC 

19.0 

s 

ANP 

24.5 

s 

SHK 

33.2 

SS 

MAT 

37.7 

SS 

TP 

SHL 

5.7 

Love' 

TP 

NDI 

20.2 

Love1 

LAH 

22.9 

Love 

EC 

SHK 

35.5 

SS 

TP 

KBL 

26.7 

S 

EC 

SHK 

35.0 

SS 

MAT 

39.7 

SS 

TP 

HKC 

35.7 

SS 

TP 

HKC 

52.0 

SS 

EC 

MAT 

35.7 

SS 

EC 

MAT 

37.3 

SS 

TP 

NIL 

15.5 

S 

QUE 

21.0 

s 

TP 

NDI 

17.9 

Love 

TP 

QUE 

25.3 

s 

TP 

NDI 

20.2 

Love1 

NIL 

22.9 

Love 

TP 

QUE 

17.2 

S 

TP 

NIL 

16  4 

S 

TP 

ANP 

48.6 

SS 

TP 

ANP 

50.0 

SS 

HKC 

46.0 

SS 

EC 

ANP 

20.8 

s 

MAT 

36.0 

SS 

EC 

MAT 

35.9 

SS 

TP 

NDI 

15.1 

Love 

EC 

MAT 

34.6 

SS 

TP 

XAN 

24.6 

s 

TP 

HKC 

510 

SS 

TP 

GYA 

27.2 

s 

ANP 

39.7 

SS 

TP 

NDI 

10.1 

PS 

NIL 

134 

p* 

TP 

LZH 

20.1 

s 

XAN 

24.0 

s 

TP 

LZH 

18.4 

s 

TP 

NDI 

20.2 

Love' 

TP 

SEO 

54.5 

SS 

HKC 

69.7 

SS 

TP 

BAG 

47.2 

SS 

TP 

ANP 

44.8 

SS 

TP 

HKC 

40.1 

SS 

ANP 

44.2 

SS 

TP 

GZH 

51.8 

SS 

‘Holt  and  Wallace,  1089;  2this  study.  TP  Tibetan  Plateau;  EC  Eastern  China. 
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Table  2 

Velocity  models  for  Tibetan  Plateau  (  TIP  ) 
and  Eastern  China  (  ECH  ) 


Depth 

TIP 

ECH 

Depth 

TIP 

ECH 

Depth 

TIP 

ECH 

(km) 

(km/sec) 

(km/sec) 

(km) 

(km/sec) 

(km/sec) 

(km) 

(km/sec) 

(km/aee) 

0 

2.550 

3.450 

150 

4.550 

4.440 

290 

4.705 

4.668 

3.75 

3.500 

3.450 

160 

4.550 

4.467 

300 

4.708 

4.675 

10 

3.500 

3.700 

170 

4.590 

4.484 

310 

4.712 

4.682 

20 

3.700 

3.800 

180 

4.625 

4.501 

320 

4.715 

4.695 

30 

3.700 

3.850 

190 

4.658 

4.518 

330 

4.718 

4.710 

35 

3.700 

4.650 

200 

4.668 

4.535 

340 

4.721 

4.720 

40 

3.800 

4.650 

210 

4.675 

4.552 

350 

4.725 

4.730 

70 

4.600 

4.600 

220 

4.680 

4.569 

360 

4.730 

4.740 

80 

4.600 

4.550 

230 

4.685 

4.586 

370 

4.740 

4.750 

90 

4.600 

4.500 

240 

4.690 

4.603 

380 

4.750 

4.760 

100 

4.600 

4.450 

250 

4.690 

4.620 

390 

4.760 

4.770 

110 

4.600 

4.400 

260 

4.693 

4.637 

400 

4.770 

4.770 

130 

4.550 

4.400 

270 

4.697 

4.654 

405 

5.014 

5.014 

140 

4.550 

4.420 

280 

4.701 

4.660 

425 

5.050 

5.050 
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Table  3  Compressional  Velocities  of  TIP 


Depth 

Thick 

Vel 

Depth 

Thick 

Vel 

(km) 

(km) 

(km/sec) 

(km) 

(km) 

(km/sec) 

0.00 

0.00 

4.530 

280.0 

10.0 

8.470 

3.75 

3.75 

4.530 

290.0 

10.0 

8.480 

20.0 

16.25 

6.160 

300.0 

10.0 

8.486 

40.0 

20.0 

6.550 

310.0 

10.0 

8.510 

70.0 

30.0 

6.740 

320.0 

10.0 

8.540 

130.0 

60.0 

8.290 

330.0 

10.0 

8.560 

170.0 

20.0 

8.200 

340.0 

10.0 

8.570 

180.0 

10.0 

8.270 

350.0 

10.0 

8.594 

100.0 

10.0 

8.340 

360.0 

10.0 

8.623 

200.0 

10.0 

8.400 

370.0 

10.0 

8.653 

210.0 

10.0 

8.420 

380.0 

10.0 

8.691 

220.0 

10.0 

8.430 

390.0 

10.0 

8.730 

230.0 

10.0 

8.440 

400.0 

10.0 

8.770 

240.0 

10.0 

8.450 

405.0 

5.0 

8.810 

260.0 

20.0 

8.460 

9.280 

270.0 

10.0 

8.462 
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Table  4.  Relocation  of  February  22,  1980  Tibet  earthquake.  The  location 
30.506°  N,  88.583°  E,  and  the  origin  time  03:02:47.2;  After  relocation:  the  loca¬ 
tion  30.658°  N,  88.649°  E,  and  the  origin  time  03:02:43.7. 


STAT 

TPK 

TPKP 

TOBS1 

TOBS2 

TCOMl 

TCOM2 

DTI 

DT2 

AAE 

51.3 

50.0 

544.1 

547.6 

544.8 

545.5 

-3.1 

iW 
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57.8 

59.3 
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734.1 

733.4 

733.8 

-2.9 

I 

AQU 

51.3 

51.3 

604.1 

607.6 

607.9 

607.7 

-4.7 

-1.0 

BAG 

17.4 

16.1 

390.2 

393.7 

392.3 

392.2 

0.9 

BUL 

36.1 

36.0 
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712.4 

712.3 

713.0 

-0.0 

CHG 

14.8 

14.9 

207.6 

211.1 

212.3 

213.4 

-1.5 

COL 

27.9 

28.3 

700.7 

704.2 

705.5 

704.6 

-5.1 

-0.6 
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37.1 

37.8 

589.9 

593.4 

594.1 
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-5.1 

-1.0 

CTA 

26.1 

26.0 

698.8 

702.3 

702.7 

702.9 

-3.9 

-0.6 

DAV 

36.4 

34 

469.2 

472.7 

471.3 

471.5 

-3.5 

-0.1 

ESK 

34.7 

34.8 

647.5 

651.0 

652.1 

651.5 

-4.9 

-0.9 

GDH 

36.6 

36.3 

709.4 

712.9 

713.5 

712.7 

-4.0 

0.3 

GRM 

28.6 

761.4 

764.9 

764.0 

764.7 

-3.5 

-0.6 

HKC 

2.8 

5 

315.6 

319.1 

319.0 
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-4.9 
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HLW 

37.9 

34 

530.7 

534.2 

529.3 

529.4 

4.8 

1ST 

29.4 

25 

522.2 

525.7 

526.4 

526.2 

-1.0 

JER 

6.2 

499.0 

502.5 

500.3 

500.4 

1.3 

KBL 

40.8 

41.1 

233.6 

237.1 

237.9 

238.0 

-0.7 

KEY 

54.9 

55.0 

547.7 

551.2 

551.4 

551.2 

-4.0 

0.4 

KTG 

46.8 

46.2 

659.6 

663.1 

662.0 

663.1 

-3.6 

0.7 

LOR 

24.3 

25.1 

637.1 

640.6 

642.3 

640.6 

-4.7 

i 

o 

oo 

MAT 

32.0 

31.3 

464.8 

468.3 
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467.9 

-3.6 

0.9 

MUN 

41.9 

40 

654.7 

658.2 
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658.7 

-2.5 

0.3 

NAI 

43.9 

44 

596.7 

600.2 

599.0 

599.8 

-3.4 

NDI 

11.6 

11.5 

144.4 

147.9 

146.6 

147.6 

-1.4 

1.0 

NIL 

55.2 
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191.5 

190.9 
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-2.8 
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NUR 

48.8 

49.5 

541.6 
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57.4 
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-2.6 
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SHL 
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12 
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-1.3 

0.8 

SNG 
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Figure  1.  Sketch  map  of  plate  tectonics  of  China,  and  the  directions  of  the  mo¬ 
tions  of  the  plates  surrounding  China.  China  occupies  the  Eurasian  plate, 
part  of  the  Indian  plate,  and  a  small  part  of  the  Philippine  Sea  plate  (  Zhang 
et.  al  1984  ). 

Figure  2.  Isopach  map  showing  crustal  thickness  (  in  kilometers,  Zhang  et.  al 
1984)  in  China.  In  eastern  China,  a  distinct  north-northeast-trending  belt 
that  contains  a  marked  difference  in  crustal  thickness  is  parallel  to  a  subduc- 
tion  zone  to  the  east  of  the  Asian  continent.  The  crustal  thickness  in  western 
China,  the  Tibetan  Plateau  is  more  than  50  km. 

Figure  3.  Fundamental  Love  wave  synthetics  for  a  variety  of  models,  source 
mechanisms,  and  depths.  The  column  on  the  left  are  the  models  used  to  gen¬ 
erate  the  synthetics  on  their  right.  The  crustal  thickness  is  60  km,  and  the 
average  velocity  of  the  crust  is  3.6  km/sec.  The  velocity  of  the  mantle  is  4.6 
km/sec,  except  for  models  g.,  4.4  km/sec,  and  h.,  4.8  km/sec.  The  distance  is 
1000  km,  and  the  source  function  is  a  trapezoid  (  1,  1,  1  )  . 

Figure  4.  A  comparison  of  Love  wave  synthetics  for  different  models.  The  column 
on  the  right  displays  the  synthetics  derived  from  the  models  on  the  left.  The 
models  are  of  a  two-layer  crust  over  half  space  mantle.  The  shear  velocities  of 
the  upper  layer  is  3.6  km /sec,  that  of  the  lower  layer  is  3.8  km/sec,  and  that 
of  the  mantle  is  4.6  km/sec.  Except  for  model  d),  the  models  have  a  31.8km 
thick  upper  crust.  The  star  is  the  source,  and  the  triangle  is  the  receiver.  A 
strike-slip  source,  with  a  Gaussian  time  history,  half  width  1.26  seconds,  is 
used.  Source  depth  is  9.8  km,  distance  is  1000  km. 

Figure  5.  Shear  velocity  models  TIP  (  Tibetan  Plateau  )  and  ECH  (  southeastern 
part  of  China  )  derived  in  this  study  compared  to  profiles  of  the  Canadian 
shield  (  SNA  )  and  the  tectonic  western  North  America  (  TNA  ). 

Figure  6.  The  ray  paths  of  Love  wave  data  used  to  derive  the  velocity  distribution 
of  the  crust  of  the  Tibetan  Plateau,  where  the  stars  are  sources,  and  the  trian- 
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gles  are  stations.  The  numbers  are  the  differential  maximum  amplitude  travel 
times  of  the  data  compared  to  the  synthetics.  The  +  sign  means  the  model  is 
faster  in  seconds.  The  number  beside  the  source  is  event  number  of  Table  1. 
The  iso-bath  relief  contour  lines  with  altitude  of  1000  meters,  2000  meters, 
3000  meters  and  4000  meters  are  given  to  define  the  physical  boundary  of  the 
Tibetan  Plateau.  The  dark  star,  event  36,  is  the  location  of  the  P ^  wave 
source. 

Figure  7.  Synthetic  comparison  of  different  models  with  the  data  from  known 
source  parameters.  In  each  group,  the  first  trace  is  data  with  station  name 
and  event  number  in  Table  1,  along  with  the  distance.  The  second,  third  and 
fourth  traces  are  the  synthetics  of  models  TIP,  CHEN,  CHUN  respectively. 

Figure  8.  The  Love  wave  data  for  paths  shown  in  Figure  6  compared  with  the 
synthetics.  Upper  traces  are  the  data,  the  lower  ones  are  synthetics  when  the 
amplitude  are  normalized,  and  shifted  by  the  amounts  indicated  by 
Tob,  —  Teyn.  On  the  left  of  each  comparison  are  the  station  names  and  the 
event  number  given  in  Table  1  followed  by  the  distance,  source  depth  and 
T  — T 

1  A  ryn  • 

Figure  9.  (  a  ):  The  ray  paths,  dotted  lines,  of  S  wave  data,  and  the  bounce  points 
of  the  SS  wave  data  used  to  study  the  upper  mantle  shear  velocity  structure 
of  the  Tibetan  Plateau.  Circles  represent  the  source-receiver  distances  greater 
than  45°,  squares  less  than  45°,  Stars  are  sources,  and  triangles  are  stations. 
The  solid  lines  from  thin  to  thick  are  contour  lines  with  1000  meters,  thin¬ 
nest,  2000  meters,  3000  meters  and  4000  meters,  thickest,  above  see  level.  The 
number  beside  the  source  is  event  number  in  Table  1.  (  b  ):  Schematic  illus¬ 
tration  of  rays  producing  S  waveforms  ( left  ),  and  SS  waveforms  (  right  ). 

Figure  10.  The  comparisons  of  the  waveforms  of  the  S-waves  for  the  paths  of  the 
Tibetan  Plateau  with  the  synthetics  using  the  model  TIP.  "  A",  "B  ”,  "  C  ", 
"  D  ",  "  E  "  and  "  F  "  are  the  branch  names  of  the  corresponding  arrivals  ( 
Figure  11  ).  On  the  left  of  each  datum-synthetic  comparison,  are  station 
names,  distances  and  event  number. 

Figure  11.  Triplication  curves  of  the  model  TIP  and  Lyon-Caen’s  model,  and 
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Lyon-Caen’s  travel  time  data  of  S  waves  for  the  Tibetan  Plateau  paths.  The 
different  symbols  indicate  the  properties  of  the  path  (  Lyon-Caen,  1986,  Fig¬ 
ure  8  ). 

Figure  12.  55  and  5-wave  data  for  the  paths  with  midpoints  in  the  Tibetan  Pla¬ 
teau  and  the  synthetics  of  the  model  TIP.  The  S-  waves  are  aligned  with  syn¬ 
thetic  S-waves.  When  the  observed  5-waveforms  are  used  as  the  source  func¬ 
tion.  On  the  left  of  each  data-synthetic  comparison,  are  station  names,  dis¬ 
tances  and  event  number. 

Figure  13.  The  ray  paths  of  the  5-  and  55- wave  data  used  to  derive  the  upper 
mantle  shear  structure  of  the  southeastern  part  of  China.  The  stars  are 
sources,  and  the  triangles  are  stations.  Circles  indicate  the  midpoints  of  the 
55  data.  Numbers  near  stars  are  event  numbers,  same  as  in  Table  1.  Numbers 
with  "  +  "  or  "  -  "  sign  on  the  rays  of  5-  waves  are  the  differential  travel  times 
of  the  data  relative  to  the  synthetics  using  the  model  ECH.  The  symbol  -f 
means  that  the  model  is  faster  than  the  data. 

Figure  14.  5  data  for  the  paths  of  southeastern  China  and  the  synthetics  using  the 
model  ECH.  Numbers  below  the  distances  are  event  numbers  as  displayed  in 
Table  1.  They  are  lined  up  with  maximum  amplitudes. 

Figure  15.  55  and  5  data  for  the  paths  with  midpoints  in  southeastern  China,  and 
the  synthetics  using  the  model  ECH.  The  5  waves  are  aligned  with  synthetic 
5  waves. 

Figure  16.  Comparison  of  the  existing  models  for  the  Tibetan  Plateau. 

Figure  17.  Comparison  of  the  data  and  synthetics  of  models  TIP,  Ltib,  Lind,  and 
Wang’s  model.  Ltib  is  Lyon-Caen’s  mode!  for  the  Tibetan  Plateau,  Lind  is 
Lyon-Caen’s  model  for  the  Indian  shield.  Letters  on  top  of  the  data  are  the 
branch  names  (  Figure  10  ). 

Figure  18.  Source  mechanisms,  and  synthetic  comparison  with  the  long-period  P 
wave  WWSSN  data  of  February  22,  1980  Tibetan  earthquake.  The  origin 
time  is  3:2:43.7,  epicenter  locates  at  30.658°  N  88.649° E,  source  depth  is  10 
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km.  The  strike,  dip  and  rake  are  160°,  60°  and  250°  respectively.  The  mo¬ 
ment  is  8.0X1024 dyn  cm  and  the  source  function  is  a  trapezoid  (  1,  1,  1  ). 
The  darker  traces  are  the  data,  and  the  lighter  ones  are  synthetics.  On  the 
left  of  each  data-synthetic  comparison,  are  the  component  of  the  data  used  ( 
"Z",  vertical  component  ),  station  name,  the  peak  amplitude  of  the  data  ,  and 
the  peak  amplitude  of  the  synthetics  in  centimeters  using  the  moment  shown 
in  the  head  line.  The  (  +  )  symbol  indicates  the  compressional,  and  little  cir¬ 
cles  indicate  dilational.  If  the  polarity  of  the  short-period  data  is  different 
from  that  of  its  corresponding  long- period  data,  we  use  that  of  the  short- 
period. 

Figure  19.  Synthetic  comparison  with  P ^  waveform  data.  The  darker  traces  are 
the  data,  and  the  lighter  ones  are  the  synthetics.  The  beginning  time  of  the 
corresponding  data  and  synthetics  are  140  and  190  seconds.  The  data  and 
synthetics  are  lined  up  with  absolute  time,  with  NDI  synthetics  shift  left  1.5 
seconds,  NIL  synthetics  shift  left  0.5  seconds.  The  numbers  beside  the  station 
names  are  event  numbers  given  in  Table  1. 
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Jurassic  or  Younger  Triassic  Late  Paleozoic  Early  Paleozoic 


Source  Depth  8  km  Source  Depth  28  km  Source  Depth  48  km 
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Figure  3 
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SECTION  2 

Broadband  Modeling  Along  a  Regional  Shield  Path, 
Harvard  Recording  of  the  Saguenay  Earthquake 
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Abstract 

The  Saguenay  earthquake,  November  25,  1088,  is  one  of  the  first  large 
shield  type  events  recorded  by  a  broadband  -  high  dynamic  range  instrument,  the 
Streckeisen  system,  installed  at  Harvard  station  (  HRV  ).  The  event  is  sufficiently 
large  to  be  well  recorded  teleseismically  and  thus  the  source  characteristics  can  be 
determined  by  independent  means  and  considered  known.  This  allows  a  detailed 
study  of  the  propagations!  effects  along  this  path,  at  an  epicentral  distance  of  625 
km,  where  the  strengths  of  the  surface  waves  can  be  compared  with  the  crustal 
body  phases.  Broadband  modeling  using  standard  analytical  techniques  and  fiat 
layered  models  works  amazingly  well  over  the  period  range  of  0.5  to  20  seconds.  A 
detailed  strategy  for  modeling  broadband  regional  phases  is  given  in  terms  of  a 
decomposition  of  the  vertical  and  radial  seismograms  into  three  segments:  P ^  ( 
containing  Pn,  j >P„,  «P.,  PmP,  P  coupled  PL  waves  )  ;  (  containing  5.,  «S„, 

SmS ,  etc.  )  ;  and  the  fundamental  Rayleigh  waves.  Information  about  the  upper 
crust  is  obtained  from  the  fundamental  Rayleigh  waves  while  crustal  thickness 
and  velocity  gradients  in  the  mantle  are  obtained  from  P^  and  S^.  This 
particular  crustal  model  has  a  thickness  of  35  km  with  a  sharp  moho  and  a 
substantial  gradient  in  the  top  20  km  of  the  mantle,  0.01  km /sec  per  km  for  both 
P  and  S  velocities.  The  mantle  velocities,  a— 8.2  and  0—4.55  km/sec  are  slower 
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than  expected  for  a  shield  environment.  Attenuation  is  not  required  for  waveform 
modeling  or  for  absolute  amplitude  estimation. 

Keywords:  Broadband  modeling,  Generalized  Ray  Method  (  GRM  ),  Mode 
Summation  Method  (  MSM  ),  Reflectivity 

I.  Introduction 

In  recent  years,  short-period  seismology  has  concentrated  on  small  events 
and  local  crustal  structure  while  long-period  seismology  has  concentrated  more  on 
inverting  source  mechanisms  of  larger  earthquakes  and  deep  earth  structure.  Long 
period  seismograms  at  regional  distances  have  not  received  much  attention  in 
either  source  studies  or  in  studies  of  the  Earth’s  shallow  structure.  One  reason  is 
that  only  small  events,  4<m<5,  remain  on  scale  on  the  standard  long-period 
WWSSN  at  these  distances.  The  source  properties  of  these  events  are  generally 
poorly  understood,  since  they  can  not  be  easily  studied  teleseismically.  Thus,  these 
events  are  not  so  useful  in  Earth  structure  investigations.  Larger  events  produce 
visible  Pn  and  PL  waves  and  these  waveforms  have  been  used  in  earthquake 
studies  (  e.  g.  Helmberger  and  Engen,  1880;  Wallace  and  Helmberger  1882  ),  but 
the  surface  waves  of  the  events  go  off-scale.  Most  digital  systems  have, 
unfortunately,  not  performed  well  at  regional  distances  for  assorted  reasons. 

With  the  installation  of  the  Streckeisen  seismometer  and  high  dynamic 
range  digitizer,  the  broadband  nature  of  regional  phases  can  now  be  appreciated  ( 
see  Figure  1  ).  The  HRV  recording  of  the  Saguenay  earthquake,  at  an  epicentral 
distance  of  625  km,  is  one  of  the  first  such  unclipped  recordings  of  a  sizable 
earthquake  occurring  in  a  predominantly  shield  region.  This  event  is  the  largest 
earthquake  to  occur  in  over  50  years  in  Eastern  North  America  and  was  recorded 
both  locally  and  teleseismically  (  e.  g.  North  et  al.,  1888;  and  Somerville  et  al., 
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1089  ).  These  studies  provide  seismic  source  parameters  for  the  event  (  strike 
—323*,  dip  «4>5*,  rake  —78*,  moment  — 5.0X1084  dyne- cm,  depth  —28  km  ),  as  well 
as  some  constraints  on  the  time  history. 

Figure  2  gives  the  three  component  long-period  WWSSN  recordings 
obtained  by  convolving  the  broadband  data  displayed  in  Figure  1  with  the 
WWSSN  instrument  response.  The  peak  to  peak  amplitude  of  the  Rayleigh  wave 
is  about  1.1  meters.  If  the  gain  of  the  instrument  is  taken  as  the  average 
WWSSN  (  2250  ),  the  Rayleigh  wave  goes  off-scale  on  conventional  records.  The 
first  40  seconds  of  these  records  is  controlled  by  the  (P—SV)  system.  This 
wavetrain  begins  with  Pm  and  evolves  into  PL  and  has  been  called  P ^  ( 
Helmberger  and  Engen,  1980  ).  The  later  arriving  phases  begin  with  Sn,  «S„,  etc 
and  go  off-scale  somewhere  near  the  direct  S  arrival.  The  beginning  portion  of  this 
group  which  arrives  before  the  fundamental  Rayleigh  wave  will  likewise  be  called 
for  convenience. 

Figure  3  displays  a  comparison  of  the  broadband  data  and  synthetics 
constructed  from  a  flat-layered  model.  A  brief  comparison  of  these  waveforms 
indicates  that  quantitative  methods  based  on  layered  models  work  quite  well,  at 
least  in  a  shield  environment.  This  paper  addresses  some  techniques  and  strategies 
which  were  used  in  deriving  the  above  model  along  with  numerous  sensitivity 
studies.  What  can  be  learned  from  regional  broadband  seismograms  is  our  main 
concern  in  this  study. 

H.  Numerical  Methods 

There  are  many  different  techniques  available  for  generating  synthetics  for 
a  flat-layered  model  (  see  Aki  and  Richards,  1980  for  example  ).  The  reflectivity 
method  and  various  modifications  essentially  perform  a  double  integration  and 
produce  accurate  results  if  the  integration  windows  are  properly  placed  with 


respect  to  the  range  of  ray  parameter  and  frequency.  This  technique  lends  itself  to 
vectorization  and  looks  more  attractive  than  the  other  methods  with  the 
availability  of  vector  computers.  Comparison  of  synthetics  generated  with  this 
approach  and  with  generalized  rays  (  Helmberger,  1983  ),  has  been  discussed  by 
Apsel  and  Luco  (  1983  )  for  simple  models  where  it  is  possible  to  sum  enough  rays 
to  construct  the  entire  wavefield.  The  advantage  of  rays  is  that  insight  can  be 
gained  into  the  timing  of  individual  pulses.  However,  in  strong  waveguides  the 
number  of  significant  multiples  become  large  and  the  method  becomes  impractical. 
Fortunately,  the  normal  mode  method  or  mode  summation  method  (  Harkrider, 
1964  )  works  well  for  modeling  surface  waves  at  about  this  time  window,  as  we 
discuss  next. 

In  this  study,  the  earthquake  will  be  represented  by  a  simple  point  source 
dislocation,  which  can  be  computed  by  constructing  a  linear  combination  of 
responses  for  pure  strike-slip,  dip-slip  and  45*  dip-slip  sources,  as  discussed  by 
Helmberger  and  Harkrider  (  1978  ).  Synthetics  of  these  three  types  of  source 
descriptions  are  displayed  in  Figures  4  and  5,  where  the  generalized  ray  results  are 
compared  with  those  of  reflectivity  and  the  normal  mode  method.  The  model  used 
is  the  same  as  that  considered  by  Apsel  and  Luco  (  1983  ),  in  their  study  of  the 
various  methods  of  generating  synthetics.  The  bottom  traces  display  Pm,  PmP 
and  direct  P,  followed  by  pPm  ,  pPmP  and  tPm  ,  etc.  A  similar  set  of  S  phases 
occurs  near  the  surface  waves,  which  appears  to  be  included  in  the  mode  solution 
and  produce  synthetics  in  good  agreement  with  the  ray  solution.  It  is  difiicult  to 
prove  theoretically  that  the  direct  S-arrivals  and  multiple  reflected  S  phases  are 
handled  properly  in  the  modal  approach  because  of  the  neflect  of  head  waves. 
Although  the  amplitudes  predicted  by  reflectivity  and  rays  differ  somewhat,  the 
waveforms  predicted  by  them  are  almost  identical.  The  differences  in  amplitudes 
are  caused  by  an  assortment  of  analytical  approximation  in  earch  individual 


method,  however,  in  this  particular  application  it  appears  to  produce  accurate 
enough  results  for  our  purposes.  The  locked  mode  approach  (  Harvey,  1081  ) 
could  be  used  to  construct  the  entire  solution  but  this  involves  considerably  more 
computing  effort.  We  will  use  all  these  methods  in  this  study,  exploiting  the 
particular  advantages  of  each  as  appropriate. 

ID.  Modeling  Strategy 

Given  a  stack  of  10  layers  with  30  parameters,  we  would  predict  a 
synthetic  match  comparable  to  that  displayed  in  Figure  3  by  letting  each 
parameter  vary  over  10  values.  In  this  section  we  discuss  a  strategy  of  finding  a 
good  fitting  model  without  wandering  through  this  large  parameter  space.  Two 
basic  elements  are  necessary  for  this  strategy  to  work,  namely,  broadband  signals 
and  the  bag  of  seismologic&l  tools  just  discussed. 

The  strategy  has  four  stages.  First,  we  model  a  long-period  version  of  the 
data  assuming  a  single  crustal  layer.  We  find  that  the  long-period  Pm  and  S„ 
depend  only  on  the  upper-mantle  velocities  and  the  average  properties  of  the 
crust.  Five  parameters  are  sufficient  for  initial  modeling,  that  is  two  velocities  in 
the  mantle,  two  average  velocities  in  the  crust,  and  crustal  thickness.  These 
parameters  are  easily  found  by  mostly  timing  and  identifying  particular  phases, 
Pm,  Sn ,  etc.  The  next  stage,  we  add  some  layers  to  the  crust  with  the  above 
constraint  to  model  the  fundamental  Rayleigh  wave. 

In  the  third  stage  we  model  the  broadband  waveforms  by  adding  some 
upper-mantle  structure,  which  controls  the  short-period  signals  riding  on  top  of 
Pu  and  5. ,  etc.  This  is  achieved  by  assuming  the  single  layered  crustal  model  and 
working  with  rays. 

In  the  fourth  stage,  we  compute  complete  synthetics  using  reflectivity  for 
the  layered  stack  and  make  some  final  adjustments  by  trial-and-error.  A  possible 
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fine-tuning  operation,  or  fifth  stage,  would  involve  a  formal  waveform  inversion. 
We  will  not  perform  such  an  inversion  here  but  we  will  show  a  number  of 
parameter  sensitivity  studies  of  the  models  near  our  preferred  model. 


a)  Long  Period  Modeling 

The  data  used  in  this  section  contain  the  WWSSN  long-period  response, 

and  will  be  referred  to  as  conventional  long-period  data.  The  beginning  portion  of 

these  records  are  displayed  at  the  top  of  the  Figure  6,  essentially  the  and 

wavetrains.  Since  a  layer-over-half  space  model  has  proved  effective  in  modeling 

the  portion  of  these  records  it  is  natural  to  examine  the  portion  assuming 

the  same  models  and  summing  generalized  rays.  Possible  models  are  give  in  Table 

1  with  corresponding  synthetics  displayed  in  Figure  6.  The  SH  phase  arriving  just 

after  140  seconds,  as  labeled  in  Figure  3,  appears  to  be  phase-shifted  relative  to 

Pn .  This  is  the  case,  as  is  seen  by  examining  the  receiver  functions  describing  the 

vertical  (  SZ  )  and  radial  components  {SR  )  of  motion  for  an  incoming  SV  signal: 

D  _  Wt 
KSZ *  -»■ 
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and  where  a  is  compressions!  velocity,  0  shear  velocity,  p  ray  parameter  (  Helm- 
berger,  1083  ).  For  p  “1  /0M ,  the  parameter  appropriate  for  the  Sn  head  wave,  we 
obtain 

Rsz^—4).2Si  and  R^r  *~0.013  . 

Thus,  we  expect  S,  and  «5„  to  be  much  stronger  on  the  vertical  component  than 
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on  the  radial  component.  Since  these  receiver  functions  are  very  sensitive  to  the 
station  site  conditions  we  would  expect  to  see  considerable  variation  in  observed 
phase  shifting  fron  station  to  station. 


TABLE  1 


Models  2  and  5  predict  better  fits  to  radial  components,  whereas  models  6 
and  7  do  better  on  the  vertical  components.  Since  is  the  strongest  on  the  vert¬ 
ical  component  we  chose  model  7  as  a  preferred  starting  model,  noting  that  phase 
shifts  will  be  expected  when  we  add  more  layers  to  the  crust. 

The  best  fitting  model  of  the  complete  long-period  seismograms,  LPM,  is 
given  in  Table  3  and  contains  five  layers  in  the  crust  yielding  the  average  crustal 
velocities  found  in  Table  1,  Model  7.  The  synthetic  comparison  of  LPM  is  given  in 
Figure  7.  These  synthetics  were  computed  by  summing  10  modes  and  assuming  a 
one-second  triangle  source  function.  These  synthetics  are  aligned  in  absolute  time 
except  the  tangential  component  which  is  shifted  4.5  second  to  the  left,  implying  a 
slightly  faster  SH  velocity  than  SV  velocity.  Several  dozen  models  were  investi¬ 
gated  but  this  particular  model  explains  the  timing  of  most  of  the  phases  the  best. 
Note  that  SH  (  SH  )  is  nearly  nodal  at  HRV  and,  therefore,  any  source  complexity 
or  errors  in  mechanism  produce  strong  changes  in  the  synthetics.  For  this  reason 
we  have  concentrated  our  efforts  in  fitting  the  more  stable  (  P—SV  )  system  where 
the  source  can  be  assumed  known. 

b)  Broadband  Modeling 

The  source  function  was  assumed  known  in  the  previous  section  but  as  we 
move  to  shorter  periods  we  must  be  more  concerned  about  source  complexity  as 
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mentioned  in  the  introduction.  We  will  discuss  time  history  sensitivity  later,  but 
in  this  section  we  have  assumed  an  asperity  type  model  consisting  of  three 
different  trapezoids  in  which  the  Sts  are  (  0.4,  0.05,  0.25  ),  (  0.2,  0.15,  0.15  ),  and 
(  0.08,  0.31,  0.23  ),  after  Somerville  et.  al  (  1080  ).  The  broadband  synthetics 
displayed  in  Figure  8  are  appropriate  for  the  long-period  model,  LPM.  The  overall 
fit  is  promising  but  lacks  the  shorter-period  signals  arriving  on  top  of  Sm  and  *5, . 
To  model  these  we  return  to  jumming  generalized  rays.  Figure  0  gives  the  com¬ 
parison  of  synthetics  of  the  preferred  model  short-period  model,  SPM,  (  Table  3, 
Figure  11  ),  to  the  waveform  data.  The  model  was  derived  by  adding  some  shal¬ 
low  mantle  structure. 

TABLE  2 

The  shear  velocity  gradient  at  the  top  of  the  mantle  has  a  strong  effect  on 
the  Sm  and  tSm  behavior  as  displayed  in  Figure  10.  The  models  used  in  construct¬ 
ing  these  synthetics  are  given  in  Table  2.  The  blank  boxes  in  Table  2  indicate 
that  the  layers  have  the  same  parameters  as  Model  209.  Each  model  is  different 
from  Model  209  in  one  parameter,  namely  the  shear  velocity  or  thickness.  The 
numbers  on  the  left  of  Figure  10  are  model  names.  We  tried  dozens  of  models 
with  a  crustal  thickness  of  40  km,  and  concluded  that  the  crust  should  be  thinner 
to  match  the  gradients  and  appearances  of  waveforms  of  S^.  If  we  use  a  faster 
crust,  Model  218,  see  Table  2,  we  can  see  the  earlier  arrivals  of  *PmP  group  ( letter 
a  on  the  radial  component  of  the  synthetics  of  model  218,  Figure  10  )  and 
•Pm Sf}m F  group  ( letter  b  ),  and  the  tSmSSmS  group  (  letter  e  )  phase,  and  the  Sw 
waveforms  predictions  are  too  early.  Model  213,  and  215  are  similar  to  Model  209 
except  that  the  thicknesses  of  the  second  layer  are  20  km  and  5  km  respectively. 


Model  213  predicts  smaller  amplitudes  of  5„  and  s Sm  phases,  and  discordant 
waveforms  of  »SmSSmS  group  phases.  Model  215  predicts  a  nice  fit  on  the  radial 
component  of  the  S ^  wave,  but  the  vertical  component  fit  is  not  so  good.  The 
thickness  of  the  second  layer  can  also  be  estimated  from  the  data.  The  shear 
velocity  of  third  layer  can  be  inferred  from  comparing  the  synthetics  of  models 
200,  and  211.  Comparing  the  synthetics  of  the  models  200  and  212,  we  can  see 
that  Model  200  is  better.  We  also  have  tried  other  velocities  for  this  layer,  but 
4.55  km  /sec  appears  to  be  the  best.  Model  200  is  the  shear  velocity  model  of  the 
SPM,  short  period  model,  while  the  compressional  velocities  are  obtained  by 
modeling  the  broadband  PH  data  described  below. 

A  similar  sensitivity  study  of  the  P^  waveform  was  applied  and  used  to 
adjust  the  P-velocity  gradient  while  holding  the  shear  velocity  structure  fixed. 
The  best  fitting  short  period  model,  SPM,  found  in  this  fashion  is  given  in  Table 
3.  Note  that  the  lower  portion  of  this  model  can  not  be  resolved  by  this  data  and 
therefore  we  have  simply  adopted  a  model  extended  at  constant  compressional 
and  shear  velocities.  The  appearance  of  Sm  synthetic  waveforms  of  Mode  summa¬ 
tion  Method  suggests  that  this  constant  velocity  layer  may  not  have  a  shear  velo¬ 
city  as  high  as  4.775  km/sec  in  a  layer  of  thickness  of  at  least  60  km. 

c)  Fine-Tuning 

In  order  to  fit  the  surface  waves  better,  we  combine  the  long-period  model 
and  the  short  period  model  and  arrive  at  the  best  model,  MPM,  see  Figure  11, 
and  Table  3.  At  this  step,  we  fix  the  gradient  at  the  top  of  the  mantle  and  make 
some  adjustments  in  the  five-layered  crustal  model  discussed  earlier.  The  normal 
mode  method  was  applied  here  since  the  fundamental  Rayleigh  waveform  is  the 
sensitive  to  the  crustal  velocity  gradient. 
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TABLE  3 


Figure  12  displays  the  long-period  data  along  with  the  synthetics  of  three 
models,  LPM,  MPM,  and  SPM.  The  second,  third  and  fourth  traces  of  each  group 
are  the  synthetics  of  model  MPM,  SPM  and  LPM  respectively.  Models  LPM  and 
MPM  predict  a  reasonable  fit  to  the  fundamental  Rayleigh  waveforms,  and 
models  SPM  and  MPM  predict  a  good  fit  to  the  early  portion  of  waveforms. 

As  a  final  check,  we  used  the  reflectivity  method  to  construct  synthetics 
for  the  model  MPM.  The  synthetic  comparisons  are  given  in  Figure  13.  The  first 
100  seconds  of  the  synthetics  fit  the  data  well.  In  this  calculation,  we  used  a  slow¬ 
ness  window  of  0.0-0.4  and  a  frequency  band  0.004-4  Hz,  assuming  an  exponential 
decay  factor  of  50  (  Qa  8800,  and  Qfi  6200  ). 

IV.  Source  Effects 

Thus  far,  we  have  assumed  that  the  earthquake  source  parameters  were 
well  determined  by  other  studies,  namely  that  the  strike,  dip  and  rake  angles  are 
roughly  323*,  65*  and  78*  respectively.  The  source  depth  was  estimated  by 
regional  P ^  and  teleseismic  depth  phases  to  be  near  28  km  (  Somerville  et  al., 
1989;  North  et  al.,  1989  ).  The  source  time  history  is  not  so  well  known  nor  are 
the  rupture  properties.  The  latter  become  important  for  regional  phases,  especially 
at  short  periods. 

In  this  section,  we  will  discuss  what  can  be  inferred  from  these  single  sta¬ 
tion  regional  broadband  records  if  we  know  the  velocity  structure  of  the  path. 
Assuming  that  our  velocity  model  is  correct;  we  now  perturb  the  above  source 
parameters  and  test  the  sensitivity  of  the  broadband  waveforms.  Because  of  the 
amplitude  difference  between  the  weak  P ^  waves  relative  to  the  surface  waves,  we 
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ha vi  displayed  these  sensitivity  results  in  two  sections,  namely  P^  and 
together  and  a  separate  scale  for  the  larger  5^  and  surface  waves. 

In  Figure  14  we  display  the  sensitivity  of  the  strike,  dip  and  rake  angles, 
while  holding  the  other  two  parameters  fixed  at  the  original  values.  Varying  the 
strike  appears  to  change  the  overall  amplitude  and  the  ratio  without  affecting  the 
ratio  of  to  S^.  The  ratio  of  Pn  to  PL  drops  somewhat  for  strikes  less  than 
323*  and  does  not  match  the  observation  as  well.  A  similar  feature  occurs  when 
the  dip  angle  falls  below  55*.  The  original  dip  looks  the  best.  Varying  the  rake 
produces  the  most  interesting  results  where  we  see  a  change  in  the  «S.  waveshape, 
becoming  sharper  for  pure  thrust  angles.  The  radial  synthetics  for  the  X=88*  looks 
the  best,  as  determined  by  overlay,  whereas  the  vertical  fits  X«=78*  are  slighter 
better  in  terms  of  Pm-S*  amplitude  ratio. 

Figure  15  displays  the  long-period  WWSSN  data  along  with  the  synthetics 
of  different  source  mechanisms.  The  source  mechanisms  do  not  have  much  effect 
on  the  surface  waveforms,  but  do  effect  the  amplitudes.  A  smaller  rake  or  a 
smaller  dip  appears  to  fit  the  Love  waveform  data  better  than  the  original  orien¬ 
tation.  Note  that  S  (  SH  )  is  near  a  node  and,  thus,  is  susceptible  to  nonplanar 
fault  effects.  Since  this  event  is  known  to  contain  at  least  two  strong  asperities  it 
is  possible  that  they  have  slightly  different  orientations  and  this  could  cause  the 
above  distortant  behavior.  One  of  the  reasons  for  concentrating  on  the  P—SV  sys¬ 
tem  is  the  insensitivity  to  source  descriptions  relative  to  the  unknown  details 
needed  to  study  SH  as  mentioned  earlier. 

Figure  16  shows  the  source  depth  effects  on  the  P and  waveforms.  The 
absolute  travel  time  of  Pm  decreases  with  the  depth  as  expected,  as  well  as  the 
time  separation  between  Pm  and  »PU  .  The  S ^  wavetrain  does  not  seem  to  change 
shape  as  much  as  one  might  have  expected.  The  source  depth  should  be  not 
difficult  to  determine  if  the  waveforms  of  the  P ^  and  waves  are  used  for  a  well 
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calibrated  path.  The  synthetics  for  a  depth  of  28  km  fit  the  data  the  best,  which 
should  be  of  no  surprise.  Clearly,  the  depth  constraint  is  the  most  important  of  all 
the  assumptions  about  source  properties  with  respect  to  modeling.  We  could  have 
easily  adjusted  the  crustal  thickness  or  velocities  to  compensate  sPn  times,  etc. 

Source  time  history  sensitivity  is  displayed  in  Figure  17  where  various  pos¬ 
sible  time  histories  are  given  on  the  right.  These  time  histories  originate  from 
different  data  sets  from  this  event,  and  are  all  idealizations  of  complicated 
phenomena.  The  simplest  ones  were  derived  from  modeling  teleseismic  short  and 
long-period  P-waves  as  discussed  earlier,  namely  (  b  )  and  (  c  ).  Case  (  b  )  was 
used  in  the  long-period  modeling  where  one  usually  uses  a  relatively  coarse  time 
sampling  and  thus  a  simple  triangle  is  sufficient.  The  more  detailed  information 
about  source  complexity  comes  from  the  strong  motion  data.  Far-field  approxima¬ 
tions  of  the  time  histories  depend  on  direction  with  two  possible  sources  (  a  )  and 
(  e  )  proposed  by  Somerville  and  Helmberger  (  1990  )  with  case  (  e  )  preferred. 
Case  (  d  )  was  used  in  the  broadband  modeling  discussed  earlier  which  is  an 
analytical  approximation  to  case  (  a  )  (  Somerville  et  al.  1989)  .  Case  (  f  ),  (  g  ) 
and  (  h  )  are  still  more  approximations.  An  overlay  comparison  of  source  function 
(  e  )  with  the  data  appears  to  be  the  best  although  there  is  not  very  much 
difference  between  the  synthetics. 

The  point  of  this  section  was  not  necessarily  to  help  constrain  this  particu¬ 
lar  earthquake  but  to  display  the  source  resolving  power  once  these  Green’s  func¬ 
tions  are  established.  This  model  and  associated  Green’s  functions  can  be  used  to 
study  other  events  roughly  along  this  path  such  as  older  historic  events  recorded 
at  Weston,  etc. 

V.  Discussions  and  Conclusions 

Our  velocity  structure  along  the  path  from  Quebec  to  Harvard  is  distinctly 
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slower  than  proposed  in  other  studies.  For  example,  the  shear  velocity  of  the  lid  of 
the  model,  4.6-4.7  km  /sec  ,  is  slower  than  that  of  either  the  ATL  model,  4.75 
km  /tec  (  Grand  and  Heimberger,  1084a  ),  or  the  SNA  model,  4.8  km  /sec  (  Grand 
and  Heimberger,  1084b  ).  In  a  tomographic  inversion  study  for  shear  velocity 
beneath  the  North  American  Plate,  Grand  (  1087  )  concluded  that  the  upper  140 
km  of  the  crust  and  mantle  of  this  region  was  very  fast  although  he  did  find  a 
major  velocity  boundary  along  this  path  at  greater  depths.  It  seems  that  the  data 
we  used  do  not  agree  with  this  fast  structure.  This  may  be  due  to  velocity  aniso¬ 
tropy.  SNA  model  isaS/f  velocity  model,  the  model  we  derived  is  a  SV  model. 
The  travel  times  of  Sm  and  «S,  waves  of  the  synthetics  of  model  MPM  match 
those  of  the  SH  data,  see  Figure  3  although  Sn  (  SH  )  is  nodal.  The  maximum 
group  velocity  of  Love  waves  of  model  MPM  is  0.8  percent  slower  than  that  of  the 
data,  suggesting  crustal  anisotropy,  although  small.  Another  explanation  may 
arise  from  epicentral  distance  uncertainty,  since  the  velocities  are  sensitive  to  the 
distance  and  origin  time.  For  example,  if  the  epicentral  distance  were  640  km 
instead  of  625  km,  the  model  would  be  2  percent  slower  than  the  SNA  model  ( 
Figure  18  and  Table  4  ).  Note  that  waveform  fits  are  comparable  to  those  found 
earlier.  Errors  in  distances  up  to  15  km  or  more  are  possible,  if  the  location  is 
determined  by  teleseismic  means.  For  instance,  the  distance  between  the  location 
given  by  USGS  and  that  given  by  HRV  is  12  km  for  this  event.  However,  in  this 
particular  case,  the  location,  based  on  close-in  aftershock  data  recorded  by  the 
Geological  Survey  of  Canada  (  North  et  al.,  1889  )  is  probably  accurate  to  5  km. 
Finally,  perhaps  this,  path  to  Harvard,  which  runs  along  an  old  suture  zone,  is 
indeed  slow  and  this  slowness  is  caused  by  a  regional  anomaly  (  Sykes,  1978  ). 
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TABLE  4 


In  summary,  we  presented  a  strategy  for  modeling  broadband  records  at 
regional  distances.  The  basic  procedure  consists  of  breaking  the  seismograms  into 
segments  -where  the  waveforms  of  each  segment  prove  sensitive  to  a  particular 
portion  of  the  waveguide.  Three  divisions  appear  to  be  the  most  useful,  namely,  a 
)  the  P*,  segment  containing  Pm  ,  pPu  ,  «P.,  PmP,  its  multiples  and  P  -coupled 
PL',  b  )  the  S 4  segment  containing  Sa  ,  «S.  ,  SmS,  its  multiples  and  S-coupled 
PL,  c  )  the  fundamental  Rayleigh  wave  segment. 

(a)  Pm  wavetrain 

The  long-period  component  of  the  P ^  wavetrain  is  controlled  by  the  cru¬ 
stal  thickness,  the  average  crustal  velocities,  and  the  compressional  velocity  of  the 
mantle.  Its  shorter  period  contributions  are  influenced  by  the  sharpness  of  the 
Moho  transition,  the  fine  structure  at  the  top  of  the  mantle,  compressional  velo¬ 
city  gradients,  etc. 

(  b  )  Stf  wavetrain 

This  segment  is  probably  the  most  difficult  to  appreciate.  It  has  been  stu¬ 
died  less,  and  produces  the  strongest  shortrperiod  arrivals,  usually  called  the  L, 
phase.  The  beginning  portion  of  this  wavetrain  consists  of  SH  ,  and  *5.,  as 
identified  in  this  study,  and  can  be  used  to  study  the  shear  velocity  properties  of 
the  top  of  the  mantle.  Further  into  this  wavetrain,  one  finds  stronger  arrivals 
associated  with  multiple-bounce  Sm  S  and  probably  S-coupled  PL  from  the  upper¬ 
most  surface  layer.  These  phases  are  included  in  the  synthetics  of  both  methods, 
the  normal  mode  method  and  reflectivity,  but  have  not  been  subdivided  into  indi¬ 
vidual  ray  groups  for  special  study.  The  longer  period  5-coupled  PL  waves 
trapped  in  the  crustal  layer  over  a  mantle  waveguide  have  been  studied 
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extensively,  and  are  responsible  for  contaminating  teleseismic  SV  body  "waves  to 
the  extent  that  they  are  rarely  used  in  source  studies.  The  laterally  varying  sur¬ 
face  waveguide  consisting  of  the  upper  few  kilometers  of  the  Earth  is  likely  to  pro¬ 
duce  similar  types  of  semi-trapped  shorter  period  P  waves.  Such  waves  are  easily 
identified  in  local  records  (  Dreger  and  Helmberger,  1090  ),  but  have  not  been 
investigated  at  these  ranges. 

(  c  )  Fundamental  Rayleigh  wavetrain 

The  fundamental  Rayleigh  waves  are  controlled  by  the  shear  velocity  struc¬ 
ture  of  the  crust  at  these  ranges  and  are  not  particularly  sensitive  to  the  mantle 
parameters.  This  statement  is  generally  true  for  the  fundamental  Loves  waves  as 
well  (  Zhao  et  al.,  1989  ).  Normal  surface  wave  modes  are  easily  generated  for 
possible  models  and  fitting  this  portion  of  the  record  proved  relatively  easy. 

In  conclusion,  this  study  demonstrates  that  broadband  regional  phases  can 
be  modeled.  This  study,  thus,  has  several  implications  with  respect  to  source 
discrimination  and  with  respect  to  exploration  of  the  lithosphere.  Earlier  studies 
have  established  the  ease  of  modeling  the  beginning  portion  of  a  long  period 
record  in  terms  of  PM ,  p  Pn ,  and  «P„ ;  i.  e.  a  layer  over  halfspace.  Apparently,  the 
Sm,  pS*i  and  «S.  window  is  also  easily  modeled  on  the  vertical  component.  The 
phases  can  be  seen  at  low  magnitudes  where  teleseismic  signals  are  no  longer 
recorded  above  the  noise.  For  example,  the  Saguenay  earthquake  produced  only  a 
few  high  quality  teleseismic  waveforms.  Thus,  S ^  and  P^  on  standard  long 
period  WWSSN  records  can  probably  be  used  to  study  events  down  to  magnitude 
4.5  since  P*i's  have  proved  useful  down  to  magnitude  5  and  5^’s  are  considerably 
stronger.  In  addition,  the  strongest  short  period  signals  are  related  in  timing  with 
the  long  period  arrivals.  One  can  use  the  latter  as  a  guide  in  constructing  broad¬ 
band  models  (  many  layers  )  and  scattering  models  for  still  shorter  periods. 

Presently,  it  has  proven  quite  difficult  to  discriminate  small  explosions 


from  small  earthquakes  as  'well  as  the  mi  Mt  ratio  does  at  larger  magnitudes  (  see 
Taylor  et  al.,  1089  ).  One  excellent  method  of  establishing  the  identification  of  an 
event  is  by  depth  and,  therefore,  by  depth  phases.  Calibrating  crustal  models  and 
the  development  of  regional  Green’s  functions  addresses  this  problem. 

Exploration  of  the  continental  lithosphere  In  interesting  regions  has  proven 
difficult.  Molnar  (  1089  )  and  his  colleagues  have  spent  a  great  deal  of  effort 
studying  the  collision  cone  between  the  Indian  and  Asian  plates.  The  absolute 
velocities  of  P  and  S  at  the  top  of  the  mantle  is  one  of  the  key  factors  in  estab¬ 
lishing  temperatures  there.  Furthermore,  velocity  gradients  can  probably  estab¬ 
lish  whether  heating  is  taking  place  from  the  top-down  or  vice  versa  (  Zhao  et  al. 
1989  ).  Thus,  broadband  modeling  of  regional  data  from  southern  platforms  of 
the  Soviet  Union  and  transition  zones  of  western  China  can  provide  some  key 
information  and  will  be  presented  In  future  reports. 
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Table  1  models  for  Pnl  waves 


Model 

n 

2 

3 

B 

5 

6 

7 

a 

6.40 

6.40 

0 

3.65 

3.65 

““"“I 

Thick 

40.0 

40.0 

40.0 

■ 

35.0 

40.0 

a 

8.10 

8.20 

8.20 

8.20 

8.20 

8.20 

8.20 

0 

4.70 

4.70 

4.65 

4.68 

Thick 

a  is  compressional  velocity,  0  is  shear  velocity. 
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Table  2.  Models  for  sensitivity  study  of  Pnl  -  Snl,  figure  10  209 


MODEL  209  211  212  213  214  215  216  218 
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Table  3.  Model  comparison,  SPM,  LPM,  and  MPM 


LPM 

SPM 

MPM 

a 

0 

Th 

a 

0 

Th 

Or 

0 

Th 

5.07 

3 .44 

3.40 

■ 

6.18 

3.57 

8.0 

6.24 

3.61 

1 

6.36 

3.68 

8.0 

6.30 

3.60 

35.0 

6.30 

3.70 

8.0 

6.52 

3.77 

8.0 

6.52 

3.77 

8.0 

6.57 

3.80 

8.0 

6.58 

3.80 

3.0 

7.00 

4.60 

10.0 

7.00 

4.60 

10.0 

8.20 

4.68 

8.10 

4.70 

10.0 

8.10 

4.70 

10.0 

8.36 

4.775 

8.20 

4.55 

00.0 

8.20 

4.55 

00.0 

8.27 

4.723 

8.27 

4.723 

35.0 

8.27 

4.723 

35.0 

8.13 

4.64 

45.0 

8.13 

4.64 

45.0 

8.13 

4.64 

45.0 

9 
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Table  4.  Model  comparison 


625  km 

640  km 

a 

0 

Th 

a 

0 

Th 

6.04 

3.49 

8.0 

6.26 

3.56 

8.0 

6.24 

3.61 

8.0 

6.45 

3.67 

8.0 

6.30 

3.70 

8.0 

6.59 

3.75 

8.0 

6.52 

3.77 

8.0 

6.71 

3.82 

8.0 

6.58 

3.80 

3.0 

6.71 

3.82 

3.0 

7.90 

4.60 

10.0 

8.10 

4.70 

10.0 

8.10 

4.70 

10.0 

8.35 

4.80 

10.0 

8.20 

4.55 

90.0 

8.40 

4.65 

90.0 
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FIGURE  CAPTIONS 


Figure  1:  The  Broadband  displacement  data  recorded  by  Harvard  station  (  HRV  ) 
of  the  25  November  1988  Saguenay  earthquake.  The  first  trace  is  vertical 
component,  where  upward  is  positive;  second  trace  is  radial  component,  posi¬ 
tive  southward;  and  the  third  trace  is  tangential  component,  positive  east¬ 
ward.  The  numbers  on  the  right  are  the  maximum  amplitudes  in  mm.  The 
instrument  response  is  flat  in  velocity  between  0.0027  Hz  to  7  Hz.  Integrating 
within  this  band  produces  the  displacement  records  displayed. 

Figure  2:  The  broadband  data  convoluted  with  the  WWSSN  long  -  period  instru¬ 
ment  response.  These  would  be  the  real  long  period  record  with  a  gain  of 
2250.  The  three  components  are  the  vertical,  radial  and  tangential.  The 
letters  on  the  left  indicate  the  positive  directions. 

Figure  3:  This  figure  displays  the  first  hundred  seconds  of  the  broadband 
waveform  data  along  with  the  synthetics  generated  by  a  flat-layered  model 
found  in  this  study.  The  main  phases  such  as  Pm ,  »Pn,  and  PL,  SH ,  and  »Sn 
are  labeled. 

Figure  4:  Comparison  of  the  radial  components  of  the  synthetics  generated  by 
Generalized  Ray  Method  (  GRM  ),  "  rays"  in  this  figure,  reflectivity,  refl,  and 
Mode  Summation  Method  (  MSM  ),  modes.  The  numbers  on  the  left  are  the 
distances  of  the  three  nearest  traces.  The  left  column  is  for  the  strike  slip 
source;  the  middle  for  dip  slip  source;  and  the  right  for  45*  dip  slip  source. 
The  source  function  is  a  0.2  0.2  0.2  trapezoid.  The  moment  is  lCPdyne-cm. 
The  source  depth  is  8.0  km.  The  numbers  on  the  right  margin  of  each  oolumn 
are  the  amplitudes  for  the  corresponding  traces  in  centimeters.  The  parame¬ 
ters  of  the  model  are:  or— 6.2km /sec  ,/9—3.5lrm/#ec,  with  a  32  km  thick  crust; 
mantle  velocities  are  or— 8.2km  /»ec  ,0^.5km  /tec .  Mode  synthetics  include  12 
modes.  The  parameters  used  for  reflectivity  synthetics  are:  decay  factor  50, 
slowness  window  (  0.0-0.5  ),  frequency  band  (  0.0025-5  )  Hz,  Qa(80000)  and 
QA 40000). 
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Figure  5:  Comparison  of  the  vertical  components  of  the  synthetics  generated  by 
GRM  method,  reflectivity  and  MSM  methods.  The  left  is  for  strike  slip 
source;  the  middle  for  dip  slip  source;  and  the  right  for  45*  dip  slip  source. 

Figure  6:  Comparison  of  the  Jong  period  P ^  and  data  and  the  synthetics  gen¬ 
erated  for  different  models.  The  upper  group  of  traces  are  vertical  com¬ 
ponents,  the  lower  are  radial  components.  The  models  used  are  given  in  Table 

1.  The  numbers  on  the  left  are  the  model  numbers.  The  source  function  is  a 
triangle  (  1,0,1  ). 

Figure  7:  Comparison  of  the  synthetics  of  the  long  period  model,  LPM,  and  the 
long  period  data.  Ten  modes  were  used  in  constructing  synthetics  for  model 
MPM.  The  tangential  component  of  the  synthetics  has  been  shifted  to  the  left 
by  4.5  seconds.  The  upper  traces  are  the  observed  waveforms.  The  observa¬ 
tions  and  synthetics  are  on  the  same  scale. 

Figure  8:  Comparison  of  the  broadband  data  and  the  synthetics  of  the  long  period 
model,  LPM.  Eight  modes  were  used  to  make  these  synthetics.  The  tangential 
component  of  the  synthetics  has  been  shifted  to  the  left  by  1.5  seconds. 

Figure  0:  Comparison  of  the  broadband  data  and  the  synthetics  of  the  final  short 
period  model,  SPM.  The  synthetics  were  generated  with  the  GRM. 

Figure  10;  Sensitivity  study  of  the  waveforms  of  and  P^  waves  to  model 
parameters.  The  numbers  on  the  left  are  model  numbers  indicated  in  Table 

2.  The  letters  a,  b  and  c  on  the  bottom  radial  component  of  the  synthetics  in¬ 
dicate  »PmPy  §PmSPmP ,  and  tSmSSmS  groups  of  phases  respectively. 

Figure  11:  The  upper  150  km  of  the  models:  long  period  model,  LPM;  short  period 
model,  SPM  and  final  model,  MPM.  The  bottom  portion  of  SPM  is  the  same 
as  that  of  MPM. 

Figure  12:  Comparison  of  the  synthetics  of  the  different  models  with  the  long 
period  data.  The  top  traces  are  observed;  the  second  is  the  synthetics  of  the 
final  model,  MPM;  the  third  is  that  of  the  short  period  model,  SPM;  and  the 
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fourth  is  that  of  the  long  period  model,  LPM.  A  1.5  second  time  shift  to  the 
left  was  applied  to  all  the  tangential  components  of  the  synthetics. 

Figure  13:  Comparison  of  the  whole  broadband  data  with  the  synthetics  of  model 
MPM.  The  synthetics  were  generated  using  the  reflectivity  method.  The 
tangential  component  of  the  synthetic  was  shifted  1.5  seconds  to  the  left. 

Figure  14:  Sensitivity  analysis  of  waveform  changes  do  to  variations  in  strike  with 
5-65*,  X— 78*  fixed,  left  column;  in  dip  with  6—323*,  X— 78*  fixed,  middle 
column;  and  in  rake  with  0—323*,  6—65*  fixed,  right  column.  The  Green’s 
functions  used  here  are  the  same  as  those  in  Figure  3  and  13. 

Figure  15:  Sensitivity  study  of  surface  waveform  data  and  amplitude  to  variations 
in  strike,  dip  and  rake.  The  MSM  technique  was  used  to  generate  these  syn¬ 
thetics.  The  vertical  components  are  given  on  the  left  and  the  tangential 
components  on  the  right  holding  the  other  parameters  fixed  as  in  Figure  14. 

Figure  16:  Source  depth  sensitivity  analysis  where  the  upper  four  traces  display 
the  vertical  components,  and  the  lower  traces  display  the  radial  components. 
The  first  trace  of  each  set  is  the  broadband  data  followed  by  synthetics  for 
various  depths. 

Figure  17:  Sensitivity  analysis  of  the  waveform  effects  do  to  variations  in  source 
history.  The  left  column  are  the  vertical  components  of  the  data  and  synthet¬ 
ics.  The  middle  two  are  radial  and  tangential  components.  The  column  on  the 
right  indicates  the  various  source  histories  used  with  the  corresponding  syn¬ 
thetics  on  the  left.  These  are  reflectivity  synthetics. 

Figure  18:  Comparison  of  the  synthetics  generated  from  model  MPM,  and  the 
model  for  640  km  (  Table  4  )  with  the  broadband  data  displayed  as  the  mid¬ 
dle  trace  of  each  set. 
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Broad-band  Modeling  of  Local  Earthquakes 
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Broad-band  Modeling  of  Local  Earthquakes 


by  Douglas  S.  Dreger  &  Donald.  V.  Helmberger 


Abstract 

Three  component  broad-band  waveforms  of  two  small  earthquakes  near 
Upland,  California  recorded  on  the  Pasadena  broad-band,  high  dynamic  range 
instrument,  were  forward  modeled  to  obtain  useful  Green's  functions.  The 
sensitivity  of  the  synthetic  seismograms  to  differing  layered  structure,  boun¬ 
dary  sharpness,  and  two-dimensional  structures  was  also  investigated.  We  as¬ 
sumed  that  the  sources  of  these  events  were  both  simple  and  known,  as  deter¬ 
mined  from  the  Caltech-USGS  array  first- motions.  A  trapesoidal  time  func¬ 
tion  was  chosen  such  that  the  width  of  the  direct  S-wave  was  well  modeled. 
The  waveforms  were  forward  modeled  using  Generalised  Rays,  Reflectivity, 
and  Finite-Difference  techniques.  In  addition,  estimates  of  moment,  fault  di¬ 
mension,  and  stress  drop  were  computed. 

The  results  of  the  modeling  exercise  indicate  that  a  simple  layer  over  a 
half-space  model  is  a  adequate  approximation  of  the  upper  crust  along  this 
profile.  In  particular,  the  waveforms  are  controlled  by  a  relatively  slow, 
surficial  layer.  The  sensitivity  analyses  indicate  that  the  details  of  boundary 
sharpness,  and  deep  crustal  structure  are  not  very  important  at  frequencies 
less  than  1  Hz,  and  at  ranges  less  than  60  km.  The  thickness  of  the  surface 
layer  was  constrained  to  be  greater  than  3  km  and  leas  than  5  km  thick.  The 
boundary  sharpness  sensitivity  study  indicated  that  the  contact  between  the 
surficial  material,  and  underlying  material  can  be  smoothed  to  2  km  without 
adversely  affecting  the  synthetic  waveform  shapes.  A  number  of  two- 
dimensional  finite-difference  calculations  were  performed,  and  it  was  found 
that  a  ridge  structure  beneath  the  recorder,  acted  as  a  lowpass  filter,  and  the 
lower  frequency  phases  were  largely  unaffected.  Other  two-dimensional  models 
with  ridges  between  the  source  to  receiver  clearly  did  not  fit  the  data  however, 
the  extended  duration  of  these  synthetic  waveforms  may  be  important  in 
modeling  more  complicated  waveforms  observed  in  data  from  events  south  of 
this  study.  Synthetic  seismograms  computed  for  the  best  fitting  model  were 
used  to  estimate  a  long  period  moment  of  6X1022  dyne-em  (ML=  4.6)  and 
1.X1022  dyne-cm  (A/^=3.7)  with  identical  triangular  source  time  durations  of 
0.3  seconds.  Assuming  the  same  fault  dimension  of  0.4  km  from  standard  scal¬ 
ing,  stress  drop  estimates  of  410  and  70  bars  are  obtained.  Generally,  we 
found  that  it  is  possible  to  model  local  waveforms  to  frequencies  of  1  Hz., 
without  a  complete  understanding  of  fine  structural  detail.  Resulting  Green’s 
functions  can  be  useful  in  studying  historic  events,  and  in  simulation  of  large 
events  from  a  given  source  region. 

Introduction 

With  the  installation  of  the  broad-band,  high  dynamic  range  Streckeisen 
instrument  at  Pasadena,  California,  it  has  become  possible  to  compare  the 
waveforms  of  different  sized  events.  In  many  cases,  events  in  the  range  of 
Mjy=3.0  to  about  M^=  5.0,  for  a  particular  source  region  have  similar  waveforms. 
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Figure  1  shows  the  location  of  the  Pasadena  station  relative  to  two  of  these  earth¬ 
quakes,  and  the  two  small  Upland  events  which  will  be  discussed  in  detail  in  this 
paper.  These  events  are  reasonably  well  located,  since  they  occur  within  the  dense 
USGS  -  Caltech  array.  In  addition,  focal  mechanisms  from  first  motions  are  easily 
obtained,  and  appear  quite  accurate  for  events  of  this  magnitude  scale.  Larger 
events  sometimes  show  complexity  because  of  the  presence  of  asperities,  on  distri¬ 
buted  fault  zones.  Furthermore,  the  frequency  of  occurrence  of  earthquakes  in 
this  size  range  is  much  higher  than  for  the  occasional  larger  earthquakes,  giving  a 
greater  availability  of  data.  In  short,  we  can  assume  the  source  is  known  (point 
source),  and  analyze  these  seismograms  in  terms  of  propagations!  effects  and 
study  the  crustal  structure.  Once  the  crustal  structure,  or  Green’s  functions  are 
obtained,  they  can  be  used  to  study  the  sources  of  sparsely  recorded  historic 
earthquakes,  or  to  scale  up  the  source  and  estimate  strong  ground  motion  for 
larger  earthquakes. 

On  traditional  instruments,  such  as  the  short-period  Wood-Anderson 
(WASP)  and  the  long-period  Wood-Anderson  (WALP),  there  is  a  rather  narrow 
range  of  earthquake  magnitudes  which  would  provide  useful  waveforms  to  model 
at  a  given  receiver  distance.  Figure  2  shows  an  example  of  two  small  local  earth¬ 
quakes,  that  occurred  near  Chino,  California  in  February,  1989  (figure  1).  These 
tangential  component  records  were  recorded  on  the  very  broad-band  channel  of 
the  Pasadena  high  dynamic  range  Streckeisen  instrument,  and  then  convolved 
with  the  instrument  response  of  a  WASP  instrument.  This  figure  shows  that 
these  two  events,  although  different  in  size  by  more  than  one  order  of  magnitude, 
have  nearly  identical  records.  This  example  demonstrates  the  deterministic  nature 
of  the  records,  indicating  that  they  could  probably  be  modeled  by  a  deterministic 
approach.  Furthermore,  note  the  maximum  amplitudes  shown  on  this  figure. 
The  larger  of  the  events  could  be  used  if  recorded  on  a  traditional  instrument, 
however  digitization  would  be  difficult.  The  small  event  could  not  be  used. 

The  purpose  of  this  paper  is  to  demonstrate  the  methodology  of  obtaining 
useful  Green’s  functions,  as  well  as,  to  show  that  quite  a  lot  can  be  learned  about 
earthquake  sources,  and  propagation  effects  using  broad-band  waveforms  at  just 
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one  site.  This  paper  is  a  two  part  study.  First,  a  crustal  model  is  determined  and 
the  modeling  approach  'is  discussed.  Secondly,  the  sensitivity  of  the  synthetic 
seismograms  with  respect  to  two-dimensional  model  perturbations  is  examined. 

Data  and  Processing 

The  Caltech  Streckeisen  instrument  is  both  broad-band,  and  high  dynamic 
range.  The  instrument  response  is  flat  in  velocity  between  0.0027  Hz  and  7  Hz. 
Integrating  within  this  band  produces  displacement  seismograms.  These  displace¬ 
ment  seismograms  can  then  be  convolved  with  the  response  of  any  instrument  to 
obtain  the  equivalent  instrument  seismograms.  In  this  study  the  displacement 
seismograms,  as  well  as,  the  WASP  and  WALP  instruments  seismograms  are  used 
in  the  forward  modeling  approach. 

Two  small  earthquakes  (Mi  =  4.6  and  3.7)  were  recorded  by  this  instru¬ 
ment  on  June  26,  1988  and  July  6,  1988,  respectively.  Both  of  these  earthquakes 
occurred  in  the  Upland,  California  area  (figure  1),  and  the  hypocenters  located  by 
the  Southern  California  Network  are  within  1  km  of  each  other.  These  events  oc¬ 
curred  at  distances  of  about  43  km  and  the  azimuth  to  the  Pasadena  station  is 
272°.  Figure  3  displays  the  three  component  displacement  data  and  the  same 
convolved  with  the  WASP  instrument.  It  is  evident  that  these  two  earthquakes 
have  very  similar  waveforms,  however  the  similarity  is  not  as  great  as  for  the 
earthquakes  presented  in  figure  2.  Since  the  events  occurred  in  nearly  the  same  lo¬ 
cation,  the  differences  in  the  waveforms  are  probably  due  to  source  complexities. 
The  most  significant  difference  in  the  waveforms  for  the  two  events  is  in  the  rela¬ 
tive  amplitudes  of  the  various  S-wave  phases  on  the  tangential  components.  The 
identification  of  these  phases,  and  an  understanding  of  the  differences  of  their 
waveforms  is  a  priority  in  forward  modeling  these  events  to  obtain  useful  Green’s 
functions. 


Methodology 

The  approach  taken  in  this  study  was  to  identify  the  phases  in  the 
waveforms  via  forward  modeling  using  Generalized  Ray  Theory  (GRT).  GRT  is 
especially  useful  for  this  purpose  because  it  allows  one  to  build  the  waveform  with 


100 


individual  rays.  Before  the  modeling  process  can  begin,  some  of  the  parameters 
must  be  constrained  by  other  means,  or  assumptions.  The  source  parameters  (lo¬ 
cation,  focal  mechanism)  are  constrained  from  first  motion  studies  utilizing  the  en¬ 
tire  Caltech-USGS  network  data,  as  mentioned  earlier.  The  mechanisms  for  the 
two  events  are  ^—125°,  £*85°,  X-»130®,  and  ^■■305°,  ^85°,  X»— 150°  (  L. 
Jones,  personal  communication).  We  assumed  a  point  source  for  both  events,  and 
a  trapezoidal  source  time  function.  A  0.15  sec.  rise-time,  0.15  sec.  duration,  and 
0.15  sec  fall-time  was  determined  by  modeling  the  width  of  the  direct  S-wave. 
The  same  source  time  function  was  used  for  both  the  main  shock  and  the  aft¬ 
ershock.  All  of  the  synthetic  seismograms  which  will  be  presented  in  this  paper 
were  constructed  with  the  focal  mechanism  of  the  mainshock.  Initial  modeling  of 
the  records  indicated  that  differences  in  the  synthetic  waveforms  were  not 
significant  for  the  two  different  mechanisms.  With  the  source  parameters  con¬ 
strained,  the  velocity  model  was  perturbed  until  good  fits  to  the  data  were  ob¬ 
tained. 


Modeling  Results 

Forward  modeling  the  tangential  component  seismograms  with  GRT  re¬ 
vealed  that  a  near  surface  layer,  with  a  relatively  low  S-wave  velocity  was  re¬ 
quired.  The  seismograms  for  both  events  posses  a  large  amplitude  phase  which 
arrives  after  the  initial  direct  S-wave.  This  phase  is  identified  as  a  near  receiver 
multiple  within  the  near  surface  layer.  Figure  4  shows  how  the  synthetics  are 
developed  with  the  addition  of  rays,  aiding  the  identification  of  the  the  multiple 
phases  (Sj  and  52).  The  phases  S0,  Sv  S2  are  all  observed  in  the  data,  especially 
on  the  tangential  records.  Note  the  ramp  like  feature  beginning  at  the  P-wave 
traveltime  and  continuing  to  the  arrival  of  the  direct  S-wave  (S0).  This  phase  ar¬ 
ises  with  the  addition  of  near-field  terms  in  the  calculations,  and  is  evident  on  the 
tangential  displacement  records  (figure  3).  Note  that  both  the  WASP  and  WALP 
(figure  5)  instruments  filter  out  the  near-field  phase  on  the  tangential  records.  Fi¬ 
nally,  the  longer  period  phase  arriving  after  5j  is  the  Love  wave.  Five  phases 
have  been  identified  on  the  tangential  component  seismograms  using  GRT.  It  is 
rather  surprising  that  such  a  simple  model  can  give  rise  to  such  complexity.  Fig¬ 
ures  5,  8,  and  7  compares  the  tangential,  radial,  and  vertical  synthetic  seismo- 
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grams  to  the  Upland  data,  respectively.  The  amplitude  given  in  these  figures  is 
the  maximum  amplitude  in  cm.  In  each  figure  the  amplitude  of  the  displacement 
synthetic  was  scaled  to  the  amplitude  of  the  data  by  multiplying  by  a  moment. 
The  same  moment  was  used  for  the  WASP  and  WALP  synthetics.  It  is  evident 
that  there  is  good  agreement  in  absolute  amplitudes  in  the  band-widths  of  these 
instruments.  The  synthetics  were  computed  using  a  reflectivity  code  vectorized 
and  programmed  by  Mallick  and  Sen,  and  the  model  LOHSl  (table  1).  The  fits  of 
the  synthetics  to  the  tangential  data  are  exceptional,  with  good  fits  to  the  radial 
and  vertical  components  for  the  first  20  seconds  of  waveform  data.  The  later  ar¬ 
riving  signals  are  difficult  to  model  but  show  a  strong  depth  dependence.  We  plan 
to  address  this  subject  later.  This  practical  code  however  does  not  take  the  near¬ 
field  problem  into  account,  and  as  a  result  the  near-field  phase  discussed  earlier  is 
not  include  in  these  synthetic  seismograms.  It  is  important  to  note  that  although 
the  model  used  contains  deeper  structure,  the  waveform  is  really  controlled  by  a 
layer  over  a  half-space,  composed  of  the  velocities  of  the  top  two  layers  of  the 
model  LOHSl. 

The  relative  amplitudes  of  S0  and  are  modeled  by  differing  source 
depths.  Recall,  that  we  began  by  constraining  the  source  depth  via  first  motions. 
It  was  found  during  the  modeling  process  that  the  relative  amplitudes  of  the  S0 
and  51  phases  are  most  sensitive  to  the  depth  of  the  source  below  the  layer.  The 
mainshock  is  best  modeled  with  a  source  depth  of  6  km,  while  the  aftershock  is 
best  modeled  with  a  source  depth  of  8  to  9  km.  In  addition,  there  is  evidence  that 
rupture  initiated  at  depth  and  propagated  upward  (Mori  and  Hartzell,  1990).  In¬ 
troducing  the  source  shallower  than  the  layer  changes  the  character  of  the  syn¬ 
thetic  waveforms  considerably.  In  short,  there  is  strong  evidence  which  favors  a 
source  shallower  than  that  obtained  from  first  motions  or  there  was  a  distributed 
source,  with  a  significant  moment  release  at  shallower  depth. 

Sensitivity  Analysis 

The  result  of  a  simple  structure  controlling  the  waveforms  is  surprising. 
Figure  1  indicates  that  there  is  considerable  geologic  complexity  along  this  profile. 
Ignoring  details  of  faulting,  and  folding,  the  profile  still  retains  a  ridge-basin-ridge 
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structure.  Specifically,  the  energy  released  during  the  earthquakes  propagated 
through  the  San  Jose  Hills,  San  Gabriel  basin,  and  the  eastern  margin  of  the  Ver- 
dugo  mountains.  The  second  part  of  this  paper  investigated  the  question  of  the 
sensitivity  of  the  synthetic  seismograms  to  perturbations  of  the  layer  over  a  half¬ 
space  structure.  Gradients,  deep  crustal  structure,  layer  thickness,  and  two  di¬ 
mensional  velocity  models  were  tested.  Three  component  synthetic  seismograms 
for  flat  layered  models,  were  computed  by  a  reflectivity  method.  Synthetic 
seismograms  for  the  two-dimensional  studies  were  computed  with  a  finite- 
difference  approach  (Helmberger  &  Vidale,  1088). 

Figure  8  shows  the  results  of  perturbing  the  layer  over  a  half-space  struc¬ 
ture  by  introducing  velocity  gradients.  All  of  the  perturbed  models  produce  syn¬ 
thetics  resembling  the  layer  over  a  half-space  result.  This  is  especially  true  for  the 
last  two  presented.  Here,  a  gradient  is  introduced  to  smear  out  the  sharpness  of 
the  contact  between  the  two  layers.  Some  of  the  effects  to  take  note  of  are  the  re¬ 
moval  of  the  second  multiple  ( S2 ),  and  the  broadening  of  the  Love  wave.  The 
second  model  has  a  linear  gradient  from  the  surface  to  the  half-space.  The  shape 
of  the  Love  wave  is  severely  altered,  as  are  the  radial  and  vertical  component  syn¬ 
thetics.  The  general  result  here  is  that  it  is  possible  to  smear  out  the  contact 
between  the  two  layers  and  still  obtain  synthetics  which  are  good  fits  to  the  data, 
however  it  is  important  that  there  is  a  well  defined  surface  layer.  It  was  found 
that  1  to  2  km  of  smearing  at  the  layer  contact  could  be  tolerated  in  the  synthet¬ 
ics  at  these  frequencies. 

Figure  0  shows  the  results  of  including  varying  deep  crustal  structure,  and 
the  effect  of  the  velocity/thickness  tradeoff  for  the  shallowest  layer.  LOHSl  is  the 
preferred  model,  discussed  earlier.  Table  1  gives  the  details  of  all  of  the  models 
used  in  the  construction  of  these  synthetics.  First,  in  comparing  the  synthetics  as¬ 
sociated  with  the  LOHSl,  SoCal,  and  LOHS2  models  one  can  see  that  there  is 
some  variation.  The  SoCal  model  is  an  average  model  of  southern  California 
which  is  used  to  routinely  invert  for  source  location  and  focal  mechanism  from 
first  motion  data.  The  synthetics  produced  by  this  model  do  not  fit  the  data  well, 
and  one  of  the  primary  reasons  is  the  thickness  and  velocity  of  the  surface  layer. 
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This  layer  is  5  km  thick,  "with  a  velocity  of  3.18  km  a”1.  The  greatest  effect  on 
the  tangential  components  is  that  the  Love  wave  is  not  as  developed.  LOHS2  has 
the  same  deep  crustal  model  as  LOHSl,  however,  it  has  a  thinner  and  slower  sur¬ 
face  layer.  The  effects  are  that  the  Love  wave  and  the  S2  multiple  are  more 
developed.  Generally,  the  thickness  of  this  layer  must  be  less  than  5  km,  and 
greater  than  or  equal  to  3  km. 

Model  Ll  is  a  modified  version  of  a  crustal  model  for  Pasadena  determined 
by  C.  Langston  (1089),  where  the  velocity  and  thickness  of  the  surface  layer  was 
changed  to  that  of  the  one-layered  model  discussed  earlier.  Ll  contains  a  low 
velocity  zone  at  depth.  Model  L2  is  the  same  model  as  Ll  except  the  low  velocity 
zone  was  removed,  and  replaced  with  a  increasing  gradient.  It  is  apparent  that 
the  synthetic  seismograms  do  not  change  very  much.  The  important  result  here  is 
that  the  details  of  deep  crustal  structure  are  not  particularly  important  in  the 
waveforms  of  local  earthquakes  (less  than  44  km),  except  at  higher  frequencies 
(greater  than  1  Hz.). 

Figure  10  shows  the  results  of  the  two  dimensional  finite  difference  calcula¬ 
tions.  We  begin  with  the  previously  discussed  layer  over  a  half-space  synthetics 
displayed  in  figure  10a.  This  figure  shows  that  the  reflectivity  and  finite-difference 
methods  correlate,  and  is  used  to  identify  the  occurrence  of  artifacts  in  the  finite- 
difTerence  synthetics.  The  finite  difference  calculation  clearly  reproduces  the  result 
of  the  reflectivity  calculations,  however  there  are  two  artifacts  present.  These  ar¬ 
tifacts  arise  from  reflections  from  the  grid  boundaries,  however  they  occur  later  in 
time  than  the  Love  wave,  and  do  not  cause  any  substantial  distortion  of  the  syn¬ 
thetic  waveforms.  Figures  10b  and  10c  show  the  effect  of  a  ridge  structure  intro¬ 
duced  beneath  the  Pasadena  station.  The  effect  is  basically  a  low-pass  filter. 
When  an  intermediate  velocity  is  introduced  into  the  ridge  (figure  10c)  much  of 
the  characteristics  of  the  layer  over  a  half-  space  is  retained,  although  there  is  a 
reduction  of  the  higher  frequency  energy.  The  WASP  instrument  synthetics 
would  be  affected.  Since  the  ridge  structure  attenuates  the  upward  propagating 
higher  frequency  energy,  initial  downward  propagating  energy,  reflected  back  to 
the  receiver  would  be  more  important  in  the  WASP  synthetics.  This  figure 
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demonstrates  that  a  ridge  structure  beneath  Pasadena  does  not  have  a  large  affect 
at  lower  frequencies.  Figures  lOd,  and  lOe  show  the  effect  of  the  addition  of  a 
ridge  between  the  source  and  the  receiver.  The  effect  is  to  completely  change  the 
character  of  the  waveform,  especially  the  Love  waves.  Surface  waves  propagate  to 
the  ridge,  propagate  through  the  ridge  as  body  waves,  and  then  develop  surface 
waves  on  the  other  side  of  the  ridge  (Vidale  &  Helmberger,  1088).  The  secondary 
surface  waves  are  typically  higher  frequency  than  those  for  the  flat  model.  In 
fact,  the  character  of  the  waveform  resembles  those  for  shallow  sources,  within  the 
surface  layer.  The  synthetics  for  this  type  structure  do  not  fit  the  Upland  data. 
Figures  lOf,  lOg,  and  lOh  explore  the  possibility  of  a  shallow  basin  within  the  sur¬ 
face  layer.  Figure  lOf  has  a  basin  with  a  1.8  kms _1  shear  wave  velocity  and  a  2.0 
gcm~3  density.  Figures  lOg  and  lOh  differ  in  that  the  velocity  of  the  basin  is 
linear  gradient  where  V  =  1.8  +  0.4Z.  The  geometry  of  the  shallow  basin  for 
figures  lOf  and  lOg  was  interpolated  from  the  depth  to  basement  results  of 
Yerkes,  (1067).  In  these  calculations  the  layer  at  4  km  depth  was  retained  because 
of  its  importance  in  producing  the  Sj  and  $2  multiple  phases  as  well  as  the  Love 
wave.  The  result  of  these  calculations  is  that  there  is  an  amplification  of  134  % 
to  273  %  of  the  maximum  amplitude,  and  an  extension  of  the  duration  of  the 
synthetic  seismograms,  compared  to  the  layer  over  a  half-space  synthetics  (figure 
10a).  Both  of  these  effects  are  the  result  of  focusing  of  reverberating  energy  in  the 
basin  which  dives  out  at  the  shallow  basin  boundary.  The  diving  energy  however 
encounters  the  deeper  layer  responsible  for  5t  and  S2  at  near  critical  angle,  and  is 
reflected  back  up  to  the  receiver.  The  differences  in  the  synthetic  waveforms  for 
figures  lOf  and  lOg  are  the  result  of  the  effect  of  the  gradient  on  the  development 
of  a  local  Love  wave  within  the  basin.  The  lack  of  the  gradient  allows  the  Love 
wave  to  develop,  extending  the  duration  of  the  waveform.  As  discussed  earlier 
this  trapped  energy  reaches  the  basin  boundary  and  leaks  out  as  body  waves,  and 
is  reflected  back  by  the  deeper  layer  responsible  for  St  and  $2.  Note  that  shallow 
sources  are  more  effective  at  generating  surface  waves  and  thus,  to  generate 
stronger  coda  via  this  type  of  scattering.  This  may  be  the  explanation  of  the  coda 
arriving  in  the  observed  data,  for  example  see  figure  7.  It  appears  that  case  (10c) 
closely  resembles  the  actual  situation  at  the  recording  station,  where  the  main 
effect  is  to  remove  the  highest  frequencies  for  the  fiat  layer  calculation. 
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Figure  11  shows  the  finite  difference  results  calculated  from  the  structure  of 
figure  10c  for  the  radial  and  vertical  components.  These  synthetics  are  compared 
with  the  reflectivity  calculations  for  the  layer  over  a  half-space.  The  synthetic 
waveforms  do  not  change  dramatically  for  this  two  dimensional  structure,  and  the 
Rayleigh  wave  seems  to  be  the  most  affected  phase,  with  a  slight  phase  shift.  To 
identify  scattered  ridge  phases  array  measurements  would  be  required.  The 
results  of  both  the  SH  and  P-SV  two-dimensional  studies  indicate  that  these 
structures  do  not  control  the  shape  of  the  synthetic  waveforms,  except  at  higher 
frequencies,  and  in  the  phase  of  the  Rayleigh  wave. 

Discussion  and  Conclusions 

We  envision  the  usefulness  of  these  Green’s  functions  in  several  ways. 
Firstly,  they  can  be  used  to  study  less  well  recorded  historic  events  where  the  only 
seismological  waveform  data  available  is  a  pair  of  long  period  horizontal  torsion 
recording  (WALP),  as  is  the  case  at  Pasadena.  This  is  particularly  valuable  for 
studying  growth  structures  in  basins  under  compression,  where  mechanisms  are 
likely  to  show  considerable  variation.  Seoondly,  they  can  be  used  to  calibrate  the 
empirical  Green’s  function  approach  in  the  simulation  of  large  events  from  small 
ones,  Hartzell  (1979)  and  Mori  and  Hartzel!  (1990).  The  latter  paper  uses  a  small 
aftershock  of  the  same  ML*=  4.6  Upland  event  considered  above  to  study  the  plane 
of  faulting  of  the  main  event  and  it  is  interesting  to  compare  results. 

The  fundamental  assumption  in  the  empirical  Green’s  function  approach  is 
that  the  waveform  distortions  caused  by  the  path  and  site  conditions  are  shared 
by  the  mainshock  and  the  aftershock.  If  the  mechanisms  are  identical,  one  can 
use  the  aftershock  P-waveforms  to  generate  the  main  event  P-waveforms  by  con¬ 
sidering  finite  rupturing  on  the  two  possible  focal  planes,  namely  (^»125° ,  5—85° , 
X»1300,  and  ^—221°,  40°,  X—8°).  They  found  that  southwest  trending  plane 

gave  the  best  results  indicating  left-lateral  motion  which  is  consistent  with  the 
overall  motion  in  the  Transverse  Ranges  frontal  fault  zone. 

Mori  and  Hartzell  obtained  a  moment  estimate  of  4.2X1022  dyne-cm  for  the 
main  event  from  the  P-wave  amplitude  following  an  expression  given  by  Boat- 
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wright  (1080).  Their  inversion  gave  a  source  area  of  one  km 2  and  a  stress  drop  of 
38  bars. 


We  obtained  a  larger  moment,  6X1022  +  2X1022  dyne- cm  by  modeling  the 
amplitudes  of  the  displacement  data  convolved  with  a  Press-Ewing  instrument 
response.  A  considerably  large  stress  drop,  based  on  the  far-field  time  duration,  t, 
was  estimated,  by  assuming  that: 
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where  a  is  the  fault  dimension  (circular  fault),  is  the  shear  velocity  and  8  is  the 
angle  between  the  normal  to  the  fault  plane  and  the  ray  path,  i>ee  Cohn  et.  al. 
(1982)  for  details.  The  trapezoidal  far-field  time  function  discussed  earlier,  can  be 
represented  by  a  triangular  time  function  (Helmberger  and  Malone,  1975).  Where 
rW).56T  1+8T2+0.58T3.  For  a  r  of  0.3  seconds,  this  expression  yields  an  estimate 

7  mo 

of  fault  dimension  a,  of  0.4  km.  Assuming,  8cr— - —  (Kanamori  and  Anderson, 

16a3 

1975),  a  stress  drop  of  410  bars  is  obtained.  This  estimate  of  stress  drop  is  consid¬ 
erably  larger  than  Mori  and  Hartzell’s  estimate,  however  different  models  of  fault 
rupture  were  used.  Recall  that  the  same  source  time  function  were  used  for  both 
the  mainshock  and  the  aftershock.  The  moment  of  the  mainshock  was  found  to 
be  about  6  times  larger  than  that  of  the  aftershock,  and  therefore  the  stress  drop 
of  the  mainshock  should  be  large  relative  to  the  aftershock.  Mori  and  Hartzell 
report  that  rupture  on  the  southwest  trending  fault  plane  (their  preferred  fault 
plane)  propagated  upward  and  to  the  southwest  (toward  Pasadena).  Considering 
unilateral  rupture  toward  Pasadena  would  increase  our  estimate  of  fault  dimen¬ 
sion,  and  would  result  in  a  lower  estimate  of  stress  drop.  In  any  case,  clearly  a 
combination  of  numerical  and  empirical  approaches  would  be  the  most  powerful, 
where  the  longer  period  properties,  namely  moment,  effective  time  history  and 
orientation  are  modeled  numerically,  and  shorter  period  effects  namely  rupture 
properties  are  modeled  empirically. 


In  conclusion,  it  is  possible  to  interpret  local  broad-band  records  with  rela¬ 
tively  simple  models.  The  broad-band  records  facilitate  the  identification  of 
different  phases  particularly  phases  with  different  frequency  contents  (i.e.  near¬ 
field,  post-critical,  and  surface-wave  phases).  The  results  of  the  forward  modeling 
indicate  that  the  waveforms  of  the  two  Upland,  California  earthquakes  are  con¬ 
trolled  by  the  relatively  slow  material  at  the  surface  for  which  a  single-layer  struc¬ 
ture  is  a  adequate  approximation.  The  forward  modeling  approach  used  here  is 
also  useful  in  studying  source  complexities.  This  is  facilitated  by  studying  two,  or 
more,  earthquakes  in  a  given  area.  The  differences  in  the  waveforms  can  be 
understood  in  terms  of  differing  focal  mechanism,  source-time  function,  or  in  dislo¬ 
cation  complexities.  In  the  case  of  the  A/x,=  4.6  Upland,  California  earthquake  it 
was  found  that  a  significant  portion  of  the  energy  was  released  at  a  shallower 
depth,  roughly  3  km  shallower  than  the  array  location.  This  is  not  surprising  in 
that  a  strong  trade-off  in  depth  vs.  origin  time  exists.  Thus,  it  appears  that  one 
of  the  first  impacts  of  the  addition  of  local  broad-band  waveform  data  will  be  on 
depth  control. 

The  sensitivity  analysis  indicates  that  the  layer  must  be  between  3  to  5  km 
thick,  and  that  the  boundary  of  the  layer  can  be  smeared  out  in  velocity  over  1  to 
2  km.  Deeper  crustal  structure  is  not  important  at  this  distance,  however  two 
dimensional  near  receiver  structures  attenuate  the  higher  frequencies  (>  1  Hz.). 
None  of  the  closed  basin  models,  nor  the  shallow  basin  models  satisfied  the  data, 
however,  the  extended  duration  produced  by  shallow  basins  may  account  for  the 
extended  duration  of  S-wave  packets  for  earthquakes  observed  to  the  south.  For¬ 
tunately,  the  relative  insensitivity  of  the  synthetics  to  details  of  boundary  sharp¬ 
ness,  layer  thickness,  and  to  some  degree  two  dimensional  near  receiver  structure, 
shows  that  useful  Green’s  functions  can  be  found  without  knowing  the  fine  struc¬ 
tural  details. 
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Figure  Captions 

Figure  1.  Map  showing  the  locations;  the  triangle  is  the  Pasadena  station,  the  cir¬ 
cle  is  the  Chino  events,  and  the  stars  denote  the  Upland  events.  Hatched  regions 
represent  areas  of  shallow  or  surficial  basement  rocks  after  Yerkes  (1965). 

Figure  2.  Tangential  components  of  the  February  16,  1988  (Ml=3.2)  and  Febru¬ 
ary  18,  1988  (ML=4.3)  Chino  earthquakes.  The  integrated  broad-band  recordings 
have  been  convolved  with  a  short-period  Wood- Anderson  instrument  response. 
Each  seismogram  is  scaled  to  maximum  amplitude. 


Figure  3.  Three  component,  displacement,  and  short-period  Wood- Anderson  data 
for  the  June  26,  1988  (Ml=4.6)  and  July  6,  1988  (ML=3.7)  Upland  events.  Each 
seismogram  is  scaled  to  maximum  amplitude. 


Figure  4.  The  addition  of  rays  to  construct  a  synthetic  waveform,  assuming  a 
moment  of  1.0X1025  dyne-cm.  Note  the  ramp  like  near  field  term  on  the  tangen¬ 
tial  component.  Each  seismogram  is  scaled  to  maximum  amplitude. 


Figure  5.  Comparison  of  tangential  component  synthetic  seismograms  with  data 
(synthetics  were  computed  with  generalized  rays).  Moments  are  1.7X1022,  and 
7.9X1022  dyne-cm  for  the  aftershock  and  mainshock,  respectively.  Each  seismo¬ 
gram  is  scaled  to  maximum  amplitude. 


Figure  6.  Comparison  of  radial  component  synthetic  seismograms  with  data  (syn¬ 
thetics  were  computed  with  reflectivity).  Moments  are  1.5X1022,  and  8.2X1022 
dyne-cm  for  the  aftershock  and  the  mainshock,  respectively.  Each  seismogram  is 
scaled  to  maximum  amplitude. 


Figure  7.  Comparison  of  vertical  component  synthetic  seismograms  with  data 
(synthetics  were  computed  with  reflectivity).  Moments  are  7.7X1 021,  and  1.9X1022 
dyne-cm  for  the  aftershock  and  the  mainshock,  respectively.  Each  seismogram  is 
scaled  to  maximum  amplitude. 


Figure  8.  Boundary  sharpness  sensitivity  assuming  a  moment  of  1.0X1025  dyne- 
cm.  Each  seismogram  is  scaled  to  maximum  amplitude. 


Figure  9.  Flat  layered  model  sensitivity  assuming  a  moment  of  1.0X1025  dyne-cm. 
Each  seismogram  Is  scaled  to  maximum  amplitude. 
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Figure  10.  SH  finite  difference  calculations  showing  ridge  structures  and  realistic 
basins  after  Yerkes,  1065.  Dotted  region  in  lOf  represents  p=2.0  g  cm-3  and 
#=1.8  km  s-1.  For  figures  lOg  and  lOh  the  dotted  region  represents  a  linear  gra¬ 
dient  in  velocity  with  the  form  v  -=  1.8  +  0.4z.  Moment  used  is  1.0X1025  dyne- 
cm.  Each  seismogram  is  scaled  to  the  maximum  amplitude  of  10a,  except  lOfgh, 
where  their  maximum  amplitudes  are  given. 


Figure  11.  Radial  and  vertical  component  synthetics  showing  the  difference 
between  a  flat  layer  over  a  half-space  and  the  model  in  figure  10c.  Moment  used 
is  1.0X1025  dyne- cm.  Each  seismogram  is  scaled  to  maximum  amplitude. 
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WAVE  FIELDS  FROM  AN  OFF-CENTER  EXPLOSION 
IN  AN  EMBEDDED  SOLID  SPHERE 

BY  LIAN-SHE  ZHAO  AND  DAVID  G.  HARKRIDER 

ABSTRACT 

This  study  investigates  the  effects  of  explosions  in  asymmetric  source  re¬ 
gions  on  the  excitation  of  seismic  body  waves.  We  give  an  analytic  formulation 
for  determining  the  wave  fields  from  an  off-center  explosion  in  an  embedded  solid 
sphere  in  an  elastic  whole-space.  As  expected,  this  geometry  generates  shear  as 
well  as  compressional  body  waves.  The  calculated  wave  fields  show  that  the  SV 
and  SH  wave  generation  is  determined  by  the  asymmetry  of  the  source  region. 
The  results  are  compared  with  the  known  analytic  solutions  of  an  explosion  in  an 
elastic  whole-space  and  at  the  center  of  an  elastic  sphere  embedded  in  the  whole- 
space.  The  radiation  patterns  at  different  periods  for  different  parameters  of  the 
media  suggest  that  the  asymmetry  of  the  source  region  has  significant  effects  on 
shorter  period,  but  has  only  minor  effects  on  long  periods.  The  long  period  P  to 
S  wave  maximum  amplitude  results  are  in  agreement  with  that  for  explosions  in 
axisymmetric  cavities. 

INTRODUCTION 

Shear  waves  are  almost  always  observed  from  underground  explosions. 
One  can  visualize  many  ways  to  convert  explosion  P  waves  into  SV  waves.  An  ob¬ 
vious  mechanism  is  the  free-surface  or  other  near  by  boundaries.  But  the  presence 
of  SH  waves  is  especially  perplexing.  The  SH  waves  are  generally  considered  to  be 
due  to  tectonic  release  from  volume  relaxation  of  the  local  prestress  surrounding 
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the  shot  point  or  as  a  triggered  earthquake  on  nearby  fault  planes  (Press  and  Ar- 
chambeau,  1962).  Other  mechanisms  for  SH  wave  excitation  have  also  been  inves¬ 
tigated;  such  as,  the  conversion  of  source  P  waves  from  explosions  in  long  narrow 
valleys  and  the  generation  of  S  waves  by  explosions  in  anisotropic  media  and  ax- 
isymmetric  cavities  (Glenn  etal,  1985,  1986,  and  Rial  and  Moran,  1986).  We  show 
in  this  paper  that  a  significant  amount  of  shear  waves  can  be  generated  by  an 
asymmetric  source  region.  We  give  the  formulation  for  the  wave  fields  of  an  off- 
center  explosive  source  in  an  elastic  homogeneous  sphere  embedded  in  an  elastic 
homogeneous  whole-space  with  different  elastic  constants  and  densities.  We  ex¬ 
pand  the  wave  fields  for  the  problem  into  spherical  harmonics  and  use  the  welded 
interface  boundary  conditions  at  the  spherical  contact  to  find  the  coefficients  of 
the  expansion.  Initially,  the  spherical  coordinate  system  is  oriented  with  the  source 
on  the  vertical  axis,  ie.  0-0,  so  that  only  SV  waves  are  generated.  Once  the  dis¬ 
placement  field  is  obtained,  the  coordinate  system  is  rotated  for  a  more  general 
source  location  and  expressions  are  then  given  to  convert  the  wave  fields  into  P, 
SV,  and  SH  wave  potentials  and  their  respective  displacement  fields. 

THEORY 

The  geometry  and  the  coordinate  system  are  given  in  Figure  1.  The  dis¬ 
placement  resulting  from  an  explosive  source  located  at  (ro0o^o)  *n  spherical  coor¬ 
dinate  system,  in  a  homogeneous,  isotropic,  unbounded  elastic  medium  may  be 
expressed  in  the  form  (Ben-Menahem  and  Singh  1981,  p384): 

u  o“*  k  ag  (w)A  o  Vh  <$2)(k  0R ) 

for  the  receiver  located  at  {r,6,<f>),  where  is  the  spectrum  of  explosive  source 
function;  V  is  gradient;  A<$2)(kaR)  is  spherical  Hankel  function  of  second  kind; 
/b0*=w/a;  a  is  the  compressional  velocity;  a2“(X+2/i)/p  in  terms  of  the  density  p 
and  the  lame’  constants;  R  is  distance  between  source  and  receiver;  and 
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0  47r(X+2^) 

where  M0  is  seismic  moment. 

We  now  express  the  spherical  Hankel  function  A^2)(kaR)  in  terms  of  r0,  r 
and  the  angle  8  between  the  position  vectors  of  the  source  and  receiver  as 

*  4s)(k  „R  )«£  (2/+l)y,(k.r0)*  /’>(  k  0r  )P  ,(cos0) 

/-o 

where  r>r0,  Ji(kar)  is  spherical  Bessel  function;  and  .P)(cos(?)  is  the  Legendra 
function.  The  replacement  of  5  by  0  is  equivalent  to  assuming  that  the  source  is 
located  on  the  z  axis,  that  is,  60— O.  With  this  substitution,  the  displacement 
expression  becomes 

u  o“«  k  h  (^)A  o  E°(2^ +1  )j  i (k  a r  0)L  07(k  Qr ) 

<-o 

for  r>r0  and  where  L<7(kar )  is  defined  below. 

If  we  bound  the  source  medium  by  a  sphere  of  radius  r—a  ,  we  can 
assume  the  displacement  field  inside  and  outside  as 

ul“u 0-  E0(2f+I)[a/N01(k^r)+b<L07(kOir)]  (la) 

i-o 

for  r  >  r0  ,  r  <a  and 

U2“E(2/+l)[c;N07(k^r)+d/L07(ka8r)]  (lb) 

/-  o 

for  r>a,  where  kp—oj/0\  is  the  shear  velocity;  and  0*—n/p.  The  boundary  con¬ 
ditions  to  be  satisfied  at  r-*a  are 

u  !«u  2  (2a) 

er’=(U|Hr'5(“2)  (2b) 

,  where  5  is  the  stress  dyadic.  The  displacement  and  stress  dyadic  can  be 
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conveniently  expressed  in  terms  of  the  Hansen  vectors,  which  are  the  eigenvector 
solutions  of  the  force-free  Navier  equation  (Ben-Menahem  and  Singh  1981, 
Chapter  2).  The  Hansen  vectors  in  terms  of  spherical  harmonics  are 

L  Mt)-—/  rtf)P  „(»,*)+ j/  /‘(«[i(l+l)]TB  „(«,*) 

N  ?(vW+l)P  „(MW  3^+i)/  /*W)[(((+  1)]7B  m(U) 

where 


P  P  ",  (costy1 '«  *-e ,  Y  m,m 

VT^B  „  -<e  ,  (cosOJe  - 

and  where  Z=kQr  ,  t?*=^r  ,  /  /+(z  )-;/(z )  ,  and  /  f  (z  )«A  [2\z )  . 
The  corresponding  stresses  are  given  by 

J_ 

et-H(L)-2^ko|F,,3(0P„,+ru(e)|((/+l)]JBm,] 


The  functions  F  /  ,  are  given  by 

i~7)f  «*  «.  )~7f  ^ 
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Substituting  the  above  into  equation  (l),  we  have  for  the  l  th  components  of  the 
displacements 


u  „-,k  *  e +l)i(kQ,r 0)[/ (f  ,)P ,« ,  +  i/ ,-(?,) ^-P ,«.] 

Vi  a »?  »?1  at/ 

— (2/+l)b ,/  (fi)P  ,e,  -(2i+l)b ,/  /(€,) -i-^-P  ,e  , 

N2  A  1  d 

u  2,-<2(+l)c ,^-f  rM P i* , -K2l+l)c,(  rM  .« . 

+(2/  +1  )d ,/  f'  «2)P  ,e ,  4<2(  +1  )d ,  /  ,-((,)  j-  jf  ,e  , 
and  the  /  th  components  of  the  spherical  surface  tractions 

e  r  -=(u  n)=*ko12g  (^)Ao(2/+l)j/(kair0)[F ,73(£,)P /(cos0)e r  +^/7i(ti)^-p  ,(cos0)e #J 
-{21  +1  )a  ,/ijk  *,[2N  2F  fr{Vi)P  /(cos  0)e  r  +F  h{Vi)jg P  /(cos0)e  ,] 
-2(2/+l)b  ,^k  ai(F,73(ei)P  /(cos0)e  r  +P/4,  (£i)~P  /(cos^)e  #] 

« r  ’=(u  21 H2/  +1  )c , /ijk  ^{2N  2F  fiOfeJP  I  (cos0)e  r  +F  -^L-P ,  (cos0)e  #] 

+2(2/  +1  )d  /  /^k  0s(P  /73(^2)P  i  (cos^)e  r  +P  ,7,(£2)  ^-P  /  (cos0)e  #] 

where  A^-«/(/+l)  and  the  prime  ,(' ),  denotes  differentiation  by  the  functions 
argument,  and  where  ,  17— k ,  £2“*at°  ,  and  17— Jt^#a 
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Using  equation  (2),  and  comparing  the  coefficients  of  P/(cos0)er  *u 
-jgPl(cos6)e,’U  h  P i (cos0)e r  *E(u  / )*e r ,  and  -~P / (cos0)e r -E(u ; )-e we  have 

YC-D  (3) 

where 

D-,k=g(W)A0y,(ko,r0)[/r'(ei),J-/r({i),2/<,k.|Fr^i).2/‘k.//,i(€.)!r  (4) 

Cl 

*Mal»  b|#c<tdl)r  (5) 

^■//+(>)i)  /,♦'(«  /r'Ka) 

^1  »?2 

dr]  V\  Ci  dr;  r/2  s2 

2/x1ks,NJf,+,(>)l)  2/1, k„/, •£(«,)  2^k^NJF,;,te)  2 ftk„/,73(f2)  (6) 

Mikj/^fai)  2/x1k0riP’^1(C1)  *t2k  fi/Trffo)  2M2ktt//7i(^2) 

We  can  now  determine  the  coefficients  (5)  of  the  expansion  (1)  by  equation  (3). 

For  inclusion  of  this  and  rotated  versions  of  this  source  into  a  vertically 
inhomogeneous  half-space,  it  is  frequently  convenient  to  express  the  source  in 
terms  of  T5  and  S  waves  by  use  of  the  Helmholtz  resolution, 

u  (7) 

where  is  P  the  potential  and  ¥  is  the  S  wave  potential.  For  the  external 
medium,  they  can  be  expressed  as 

■*'-T-,E(2i+l)'i,*l2,(k„r)P,(coS»)  (8) 

Ko  1—0 

*—  i~  £  (2/  +1  )c  /  A  /2)(k  ^r )  -^r-P  , (cos0)e  d  (0) 

kfi  i-o  dt/ 

Generally,  if  the  source  is  located  at  (ro,0o,^o),  ^  anc*  4  *n  equation  (8)  and  (9) 
should  be  changed  to  6  and  7,  with  cos^cos0cos0o4sin$sin#ooos(<£— <f>0)  .  The  unit 
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vectors  at  the  source  are  now  eroXer/sin5,  and  ero.  The  range  of  S  is  0  to 
360°,  and  6  is  defined  as  the  angle  measured  clockwise  from  the  position  vector  of 
the  source  to  the  position  vector  of  the  receiver.  The  potentials  can  be  written  as 
(Appendix  A): 

—  £  £  l-^«m  cosm  4>+Bnm  sinm  0]P ™(cos6)h  ^(k  ftr )  (10) 

*  a  *  “0  m  “0 

- n-'sT' «(k  ,r )  (11) 

*p  n-0  m— 0 

where  ,  (j«*l,2,3)  are  the  Cartesian  components  of  the  S  wave  potential  ^  ,  and 
A„  rn  —  (2n  em  P*m(cos0o)cosm  0O 

— {2n  +1)j~T— os^0)sinm  0O 
(n+m  j! 

Q(m}=— y-(2n  jj  ern  cn  [-^-Pncos^Km^cosm  0O 

— mP”(cos0o)coi  0ocos0osinm  0O] 

D !,}}=— ~{2n cw c„  [-^-F^cos^in^osinm 0O 
+mP™(cos90)cot  #ocos0ocosm  0O] 

C„{m)l=-y-(2n  ern  c»  [^Pnm(cos0o)cos0ocosm  0O 

+mP™(cos9Q)cot  ^osin^osinm  0q] 

y-{2n  cm  c*  I'5irP"m(cos^)cos^sinm 

—mP™(cos90)cot  0osin0ocos  m  0^] 

— -y-(2n  c„  P„m(cos0o)sinm  0O 
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Dim,am- -j-(2n  c,  Pficos90)cosm  <t>0 

where  £m*=l  if  m  -=0  ,  em  *2  if  m  ^0  .  is  the  SH  wave  potential. 

From  the  equations  (11),  we  can  see  that  SH  waves  will  be  generated  if 
C$  and  Dj3)  are  not  zero  in  Cartesian  coordinate  system.  Harkrider  et  al.  (1991) 
gives  a  detailed  derivation  of  the  surface  waves  in  layered  media  due  to  a  tran¬ 
sparent  source  expressed  in  the  general  form  (10)  and  (11). 

When  the  explosion  source  is  on  the  z  axis,  as  in  equation  (9),  there  are 
only  SV  waves  because  of  symmetry.  As  one  moves  the  source  off  the  z  axis,  the 
and  D n^)  are  no  longer  zero  and  we  have  SH  waves  as  one  would  expect. 

Even  though  equations  (10)  and  (11)  are  convenient  for  calculating  surface 
waves,  they  are  awkward  for  evaluating  the  whole  space  wave  fields.  We  intro¬ 
duce  an  alternative  approach.  If  the  source  is  located  at  (r0,0,0),  the  displacement 
fields  and  the  potentials  are  given  by  equations  (1),  (8)  and  (9).  If  the  source  is 
located  at  (ro,0o,$o),  as  discussed  above,  the  coordinate  system,  with  ero  as  north 
pole,  has  unit  base  vectors  (er,e4,e.j)  ,  and  coordinates  (r ,6,' y)  ,  where  the 
definitions  of  8  and  e7  are  given  above.  In  this  coordinate  system,  the  displace¬ 
ment  fields  and  5  potential  can  be  written  as: 

u  “(up  -busy  )e  r  +(up  4-ugv  )e  (,  (12) 

and 

*-V*r  (13) 

where 

e  t"®  r0*e  r  /s'm5 


e  l-e  7Xe  r 
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and 


“^er(er'e'^r^ 

Thus,  the  displacement  fields  in  the  coordinate  system  (r  ,6,4>)  can  be  written  as: 

u*urer+u,e|+u*e^)  (1 


tl r  “Up  +u£v, 


u  (up +u^r  )(sin0ocos0coe(d> — $o)— cos^sm#)/ sin  6, 


u  ^*=(up^ +ugy  )sin^osin(0— $0)/sin<5. 
And  the  S  potential  can  be  written  as: 

e  r  +4fe  t+ilfe  + 


V- o, 

xl&i  n0osin(<£— <f>0)  /  sin  6, 


itfs\n00cosdcos(<t>—  <f>0)— oos9tfs\n0)/s\n8. 

When  <5=0orl80°,  equations  (14)  and  (15)  no  longer  hold.  This  is  when  the 
source  and  receiver  are  on  a  line  through  the  sphere  origin  and  then  we  should  use 
equations  (8)  and  (9). 

Before  going  on  to  numerical  examples,  there  is  one  further  point  to  dis¬ 
cuss.  The  observed  seismic  moment  is 


W 

Mot,  - - - - M0 

PlQl+ir(P202~PlPl) 
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Where  M0  is  the  input  moment.  This  equation  can  be  obtained  by  using  the  equa¬ 
tion  (B4)  of  Appendix  B  and  letting  u;  go  to  zero.  Since  nuclear  explosion  sources 
are  frequently  characterized  by  A  0  or  the  more  familiar  notation,  'I'(oo),  which  is 
the  asymptotic  value  of  the  reduced  displacement  potential  for  large  time,  and 

47rpa2'I'(oo), 


we  can  write  equation  (16)  as 


1+ 


'I'(oo) 


(17) 


where  the  elastic  rigidity  Thus  the  difference  between  observed  and  shot 

point  ty(oo)  is  independent  of  density  contrast  and  in  source  media  with  similar 
Poisson  ratios  is  only  a  function  of  the  rigidity  ratio  of  the  material  in  which  the 
'I'(oo)  is  measured  end  the  shot  point  material.  The  smaller  the  source  region  rigi¬ 
dity,  the  smaller  the  observed  'fr(oo)  while  large  source  rigidities  cause 
amplification.  Since  the  minimum  realistic  ratio  of  P  to  S  wave  ratios  is  V2, 
which  corresponds  to  a  Poisson  ratio  of  zero,  the  maximum  possible  amplification 
is  a  factor  of  3. 


RESULTS 

In  this  section,  we  present  some  numerical  results.  For  all  of  the  synthetics, 
we  assume  that  the  moment  time  history  is  a  step  function,  ie. 


and 


y(w)-i/iu;, 


MQ—AirpiCtxYXCPdync  —cm , 

where  p  is  in  g /cm3,  and  the  velocities,  a  and  0,  are  in  bn  /sec,  and  distances  in 
km .  Thus  'I'(oo)  is  kept  the  same  for  the  different  spherical  source  media  models. 
The  distances  are  in  km  and  the  displacement  is  in  cm.  The  spectral  band-width 
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is  from  102.4  seconds  in  period  to  10  Hz.  As  our  standard  model,  we  assume  the 
elastic  parameters  and  dimensions  given  in  Table  1.  The  receiver  is  located  10  km 
from  the  center  of  the  elastic  sphere  unless  otherwise  indicated.  The  material 
model  for  the  whole-space  outside  of  the  sphere  is  reasonable  for  sedimentary 
regions,  see  for  example  Zhao  et  al.  (1991).  When  investigating  the  effect  of 
different  elastic  properties  of  the  sphere,  Nafe’s  velocity- density  relation  (Ludwig 
et  al.  1970)  is  used  to  obtain  densities  and  S  wave  velocities  for  a  given  P  velocity. 
We  used  20  terms  to  sum  up  the  wave  fields,  see  Equation  (1).  The  error  induced 
from  truncation  of  order  m,  is  approximately  the  order  m6(r0/r)m. 

Figure  2  gives  a  comparison  of  displacement  fields  calculated  by  our  infinite 
series  formulation,  ie.  inversion  of  equation  (5),  and  by  the  numerical  evaluation 
of  two  analytic  closed  form  solutions.  The  upper  pair  of  synthetics  show  the 
whole  space  results,  while  the  lower  pair  show  the  results  for  a  source  located  at 
the  center  of  the  solid  sphere  in  a  whole-space.  The  formulation  for  the  wave  field 
outside  of  the  solid  sphere  due  to  a  centered  explosion  is  given  in  Appendix  B. 
The  solutions,  both  closed  and  series,  are  calculated  in  the  frequency  domain,  and 
then  transformed  to  the  time  domain,  so  that  each  pair  of  time  histories  show 
almost  identical  precursor  artifacts,  which  are  caused  by  the  numerical 
transforms.  The  agreement  within  pairs  is  excellent. 

In  the  following  figures,  we  will  refer  to  P  and  S  wave  arrivals  as  radial 
and  azimuthal.  These  are  rip,  rip,  rig,  and  rig  respectively  in  the  notation  of  the 
previous  section.  The  P  waves  are  the  contribution  to  the  displacement  field  from 
the  P  potential  , begirt  at  the  P  wave  arrival  time,  and  physically  correspond  to  P 
or  compressional  waves  arriving  at  the  receiver  which  have  undergone  multiple 
refections  and  conversions  between  P  and  S  in  the  spherical  source  region.  An 
analogous  statment  applies  to  the  S  waves  with  P  replaced  by  S  and  compres¬ 
sional  by  shear.  At  a  range,  r,  of  10  km  for  the  standard  model  ,  the  near-field  P 
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arrival  t ip,  which  is  the  displacement  component  perpendicular  to  the  propagation 
direction,  and  the  near-field  S  arrival  ■which  is  the  component  in  the  radial  or 
propagation  direction,  are  about  10%  of  the  far-field  arrivals  tip  and  «/.  For  a 
dynamic  point  source,  one  would  expect  far-field  dispacements  to  fall  off  as  1/r 
and  the  near-field  to  attenuate  as  1/r2.  This  was  investigated  by  generating  wave 
fields  at  ranges  from  5  to  40  km  in  5  km  intervals  for  the  standard  model.  The 
wave  form  of  each  of  these  displacement  components  was  found  to  be  independent 
of  range  and  by  10  km,  the  displacements  fields  indeed  attenuated  as  predicted. 

In  most  of  the  following  examples,  6  is  set  to  90°,  ie.  the  source  and  receiver 
are  both  in  the  horizontal  plane.  This  was  done  to  check  our  numerical  algo¬ 
rithms  for  the  transformations  given  in  equation  (14),  which  are  used  to  obtain 
SV  and  SH  displacements  from  the  initial  SV  only  shear  displacement  field. 

Figure  3  shows  displacements  at  a  receiver  located  at  (r  ,#,$)“ (10,90°, 90c), 
with  the  source  at  (0.45, 90°, 180°).  At  these  coordinates,  the  source  and  receiver 
position  vectors  are  90°  to  each  other  with  respect  to  the  sphere  origin.  We  see 
that,  strong  shear  waves,  as  evidenced  by  the  azimuthal  component  of  S  waves, 
are  generated  at  receiver  azimuths  90°  from  the  symmetry  axis  of  the  source  and 
the  center  of  the  sphere.  Abo  the  time  difference  between  the  two  biggest  arrivals 
in  both  the  P  and  5  waveforms,  is  about  2.5  seconds,  which  is  the  time  needed 
for  an  S  phase  to  travel  through  the  sphere.  This  interval  can  be  detected 
between  many  arrivals.  Thus  the  bigger  arrivals  in  the  seismograms,  except  for 
direct  P  and  S  generated  at  the  sphere  boundary  can  be  interpreted  as  phases, 
which  have  traveled  through  the  sphere  at  least  once  as  S.  The  S  wave  displace¬ 
ment  is  about  a  factor  of  10  greater  than  the  P  wave  displacement.  This  is  some¬ 
what  misleading,  since  as  we  shall  see  in  later  examples,  this  is  the  azimuth  of 
minimum  P  displacement  and  of  maximum  long  period  S. 

The  change  in  displacement  fields  as  a  function  of  source  asymmetry  is 
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displayed  in  Figure  4  .  The  model  and  the  geometry  are  the  same  as  that  used  for 
Figure  3,  except  that  the  source  radius,  r0  is  varied  from  0  to  0.45  km,  the  loca¬ 
tion  in  Figure  3,  in  0.05  km  increments.  From  this  figure,  we  see  that  only  P 
waves  are  generated  when  the  source  is  at  the  center  of  the  sphere,  r0  equal  to  0. 
At  this  azimuth,  as  r0  increases,  the  P  wave  maximum  amplitude  is  reduced  by 
about  a  factor  of  10,  while  the  S  wave  amplitude  is  increased  to  almost  50% 
greater  than  the  P  waves  due  to  the  centered  source.  The  bottom  traces,  r0— 0.45 
km  ,  are  the  same  as  those  in  Figure  3  except  on  a  common  scale. 

Figure  5  illustrates  the  effects  on  the  wave  fields  as  the  the  velocity-density 
structure  of  the  spherical  source  region  is  varied.  The  exterior  medium  properties 
are  same  for  all  cases  and  are  given  in  Table  1.  The  first  row  of  traces  are  the 
same  as  those  in  Figure  3.  For  the  bottom  row,  the  sphere  and  surrounding  space 
have  identical  properties  and  thus  these  traces  are  the  results  for  an  explosion  in  a 
whole  space.  The  source-receiver  geometry  is  the  same  as  Figure  3.  We  see  that, 
the  greater  the  elastic  contrast  between  the  spherical  source  region  and  exterior 
medium,  the  greater  the  S  waves  generated.  This  is  evidenced  by  the  relative  exci¬ 
tation  between  P  and  S  on  any  given  row.  In  making  amplitude  comparisons 
between  differing  source  media  in  the  figure,  one  must  remember  that  the  moment 
increases  as  one  moves  down  the  row  with  increasing  velocities  and  density  while 
'J'(oo)  remains  unchanged. 

Figure  6  shows  the  wave  field  variations  as  the  angle,  b,  between  the  source 
and  receiver  position  vectors  changes.  When  6—0,  the  source,  the  receiver,  and  the 
center  of  the  sphere  are  in  a  line.  At  this  angle,  no  P  or  S  6  components  are  gen¬ 
erated.  The  near-field  S  r  components  are  P  to  5  conversions  at  the  boundary 
which  are  focused  back  on  the  axis.  The  elastic  sphere  acts  as  a  contact  lens  mag¬ 
nifying  the  the  wave  fields,  like  the  Earth  does  at  the  antipode  of  an  earthquake 
(Rial,  1980).  When  b  increases,  the  P  amplitude  decreases,  and  5  amplitude 
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increases. 

The  bottom  traces  are  again  those  of  Figure  3.  Comparing  the  bottom  row  with 
the  top  row  and  then  with  the  top  row  of  Figure  4,  we  see  that  the  maximum  S 
wave  amplitude  is  almost  3  times  the  maximum  P  wave  amplitude  for  this  the 
maximum  degree  of  asymmetry  investigated  and  is  more  than  50%  greater  than 
the  P  wave  amplitude  for  the  center  explosion.  Thus  from  Figures  4,5  and  6,  we 
see  that  the  degree  of  asymmetry  of  the  source  region  and  the  contrast  between 
the  media  determine  the  intensity  of  the  5  wave  fields. 

In  Figure  7,  we  show  the  results  for  different  sizes  of  the  elastic  sphere.  The 
source  radius  is  kept  at  the  same  fraction  of  the  sphere’s  radius  as  in  previous 
cases,  ie.  r0  is  0.9a.  The  fourth  row  in  the  figure  is  the  same  as  the  last  row  in 
Figure  6  since  the  radii  and  source-receiver  angles  are  the  same,  ie.  0.5  km  and 
90°.  Most  of  the  observed  effects  can  be  explained  by  noting  that  increasing  the 
radius  under  these  conditions  is  similar  to  decreasing  the  velocity.  Thus  we  expect 
time  differences  between  arrivals  to  increase  and  a  decrease  in  frequency  content  as 
the  radius  increases. 

Figure  8  displays  the  radial  components  of  the  P  waves  on  the  symmetry 
axis  for  different  source  radii  with  6-0,  ie.  source  and  receiver  on  the  same  side  of 
sphere  center  and  180°,  ie.  opposite  sides.  The  distance  used  is  10  km  from  the 
source  while  in  previous  figures  the  receiver  was  10  km  from  the  sphere  center.  If 
the  source  is  located  at  the  center  of  the  solid  sphere,  the  wave  fields  for  6-0,  and 
180°  are  of  course  identical.  The  numerical  evaluation  of  the  analytical  solution 
was  shown  in  Figure  2.  As  the  source  is  moved  toward  the  boundary  of  the 
sphere,  the  early  arrivals  separate  and  diminish  in  amplitude  for  6-0.  For  6—180° 
,  they  converge  and  grow  in  amplitude  as  the  source  passes  through  an  effective 
focus  in  the  sphere  and  then  decrease  in  amplitude  as  the  source  is  moved  closer 


139 


to  the  boundary. 

In  Figure  9,  we  keep  the  dimensions  and  source-receiver  angles  the  same  as 
Figure  8  but  vary  the  sphere  elastic  properties  using  the  same  relationship 
between  velocities  and  densities  as  before.  The  source  radius  is  0.45.  The  first  row 
is  the  radial  P  wave  field  for  velocities  and  density,  which  are  the  same  as  the 
previous  figure.  They  can  be  considered  an  additional  row  to  the  previous  figure 
with  the  source  radius  increased  from  0.4  to  0.45  km.  All  traces  on  this  figure  are 
on  a  common  scale  in  order  to  demonstrate  the  increase  in  amplitude  of  the  early 
arrivals  at  6—0°  and  the  reduction  in  amplitude  of  the  later  arriving  S  coupled 
internal  reflections  at  both  receiver  positions  as  the  velocities  and  density  of  the 
sphere  are  increased  to  that  of  the  whole-space.  The  whole-space  result  is  only 
slightly  larger  than  the  6-0°  trace  on  the  bottom  and  of  course  smaller  than  the 
6*180°  bottom  trace.  The  effect  of  decreasing  the  impedance  contrast  and  of 
reducing  the  travel  time  for  reflections  in  the  sphere  reduces  the  partitioning  of 
direct  P  into  later  arrivals.  The  increase  in  moment  while  keeping  Moo)  constant 
as  we  go  down  the  rows  seems  to  have  only  a  minor  effect  on  the  direct  P  ampli¬ 
tude.  It  is  notable  that  the  largest  amplitude  of  direct  P  is  in  the  middle  rows 
where  the.  elastic  contrast  is  intermediate  and  in  the  receiver  direction  where  the 
path  is  longest  in  the  sphere,  ie.  £=180°.  This  is  due  to  the  source  being  near  the 
focus  which  changes  position  relative  to  source  point  as  we  change  velocities.  This 
is  similar  to  the  previous  figure  where  we  moved  the  source  relative  to  the  fixed 
focus. 

With  the  exception  of  Figure  6,  we  have  investigated  the  time  histories  at 
either  the  long  period  maxima  for  the  S  waves,  ie  6*90°  or  at  the  P  wave  max¬ 
ima,  ie.  6«=0°  or  180.  Even  in  Figure  6,  6  ranged  only  from  0  to  90°.  In  the  next 
three  figures,  we  show  the  radiation  pattern  for  P  and  5  at  periods  of  1  and  20 
seconds  as  a  function  of  source  radius,  elastic  parameters  and  radius  of  the  sphere. 
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The  periods  were  chosen  for  their  obvious  significance  to  seismic  magnitudes.  The 
asymmetry  of  the  source  region  is  seen  to  have  significant  effects  on  short  period 
radiation  patterns,  but  to  have  only  minor  effects  on  the  long  period. 

Figures  10,  11  and  12  show  radiation  patterns  for  different  configurations. 
"P"  and  "S"  in  these  figures  denote  the  radial  component  of  P  waves  and  the  S 
component  of  S  waves  respectively.  The  distance  used  for  these  calculations  is 
1000  km  to  void  the  near-field  contributions.  The  source  is  on  the  horizontal  line 
for  all  these  figures.  The  direction  of  6  'is  indicated  on  the  diagram  of  the  sphere  at 
the  top  of  Figure  10  with  a  star  marking  the  source.  The  sphere  center  is  shown 
as  a  "+"  in  all  the  radiation  patterns.  The  dashed  circles  with  the  same  centers  as 
the  P  radiation  patterns  are  the  P  wave  radiation  pattern  for  the  source  located 
at  the  center  for  the  particular  model  and  demonstrate  the  amplitude  effect  due  to 
source  asymmetry.  The  dashed  line  part  of  the  S  wave  radiation  pattern  indi¬ 
cates  negative  amplitude  compared  to  the  solid  line.  The  number  at  the  right 
hand  side  of  the  S  wave  radiation  pattern  is  the  multiplyer  of  the  S  wave  pattern 
used  to  plot  on  the  figure,  ie.  the  larger  the  number  the  smaller  the  actual  pattern 
compared  to  other  S  wave  patterns  on  the  figure  and  the  P  wave  pattern  next  to 
it. 

Figure  10,  using  the  standard  model,  shows  the  radiation  patterns  for 
different  source  locations,  r0.  Figure  12  uses  the  same  elastic  model  but  different 
radii  of  the  solid  sphere.  In  Figure  11,  the  sphere’s  elastic  parameters  are  varied  in 
the  inverse  order  to  Figure  9.  Similar  effects  are  seen  in  all  three  numerical  experi¬ 
ments.  As  the  source  radius  is  increased,  the  elastic  parameters  of  the  sphere  are 
decreased,  or  the  sphere  radius  is  increased,  the  radiation  patterns  become  more 
asymmetric.  In  other  words,  the  P  wave  radiation  becomes  less  spherical  and  the 
20  second  S  wave  pattern,  which  'is  essentially  the  radiation  pattern  for  a  point 
force  pointed  away  from  the  source,  becomes  larger.  From  these  figures,  we  also 
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see  that  the  intensity  of  the  P  waves  radiated  by  an  off-center  explosion  is  greater 
than  that  by  the  corresponding  center  explosion  with  the  same  ♦(oo)  except  in 
Figure  12  where  increasing  the  source  and  sphere  radius  beyond  the  dimensions  of 
the  standard  model,  the  opposite  is  true.  Also  by  comparison  wit'inu  each  row,  we 
can  see  that  S  wave  generation  is  much  smaller  and  that  the  P  wave  fields  are 
closer  to  the  center  explosion  excitation  for  longer  periods  ther  *•  ><*  shorter  periods. 

An  interesting  feature  is  the  shape  of  the  radiation  patterns  at  short 
periods.  For  the  P  wave  in  Figures  10  and  11,  we  see  that  the  radiation  pattern  is 
similar  to  a  dipole.  The  radiation  patterns  are  reminiscent  of  those  of  far-field  P 
and  S  waves  scattered  from  localized  inhomogeneities  when  plane  P  waves  are 
incident  {Aki  and  Richards,  1981,  pg.  732).  This  should  not  be  surprising  since  in 
our  case  the  incident  wave  is  a  spherical  P  wave  interior  to  the  scatterer.  In  com¬ 
paring  the  Aki  and  Richards  figure  with  ours,  one  must  remember  that  our  P 
wave  radiation  patterns  include  the  transmitted  spherical  wave  perturbed  by  the 
scattering  or  asymmetric  contribution  where  theirs  is  only  the  scattered  wave. 

Another  interesting  feature  in  all  three  figures  is  that  the  maximum  long 
period  20  second  P  wave  amplitudes  are  at  least  twice  («  1.8)  as  large  as  the 
corresponding  maximum  5  wave  amplitudes  for  all  cases  calculated.  This  is  in 
agreement  with  the  empirical  formula  for  explosions  in  axisymmetric  cavities  of 
Rial  and  Moran  (1986).  But  because  their  sources  were  in  the  symmetry  center  of 
the  cavity,  the  resulting  5  wave  radiation  was  more  like  a  dipole  than  a  point 
force  at  low  frequencies. 

In  Figure  13,"  we  investigate  the  error  introduced  if  we  estimate  the 
moment  by  using  a  finite  time  length  which  may  include  only  the  initial  few 
pulses.  The  synthetics  in  Figure  13  are  calculated  for  our  standard  velocity  and 
density  mode!  keeping  only  1/r  terms,  ie.  the  far-field  approximation.  The  dis¬ 
tance  used  is  10  km  from  the  center.  The  top  trace  is  the  response  of  a  centered 


explosion;  the  lower  one  is  the  response  of  an  off-center  explosion.  The  source  is 
0.45  km  from  the  center  of  the  sphere  and  the  angle, 8,  between  the  source  and 
receiver  is  45°.  From  the  long  period  P  wave  radiation  patterns,  we  know  that  the 
angle  is  irrelevant  as  long  as  we  are  interested  in  moment.  The  moment  is  calcu¬ 
lated  from  the  area  using  the  zero  frequency  spectra  and  corrected  from  observed 
to  actual  moment  by  equation  (16).  From  the  figure,  the  longer  the  time  interval 
we  use,  the  better  the  estimate  we  get.  The  first  few  pulses  hold  more  than  80  per¬ 
cent  of  the  moment  in  these  two  examples.  Even  though  the  first  pulse  holds 
about  90  percent  of  the  moment,  in  a  real  situation  it  would  be  difficult  to 
separate  the  effect  of  local  reverberation.  In  this  case,  if  we  take  a  long  enough 
window,  the  reverberations  finally  die  down  and  we  get  a  good  estimate. 

Figure  14  shows  Mobt/M0  as  a  function,  equation  (16),  of  the  sphere 
compressional  velocity, alt  for  different  values  of  the  whole-space  compression  velo¬ 
cities,  q2,  assuming  Nafe’s  velocity-density  relation.  For  exterior  velocities  greater 
than  those  of  the  sphere,  the  observed  moment  Mobt  is  greater  than  the  input  M0. 
And  for  velocities  less,  the  opposite  is  true. 

CONCLUSIONS 

In  this  study,  we  gave  the  analytical  formulae  for  the  wave  fields  from  an 
off-center  explosion  in  an  embedded  sphere.  Numerical  calculations  were  compared 
with  the  known  analytic  solutions.  We  also  calculated  P  and  5  wave  fields  for 
different  r0,  the  source  distance  to  the  center  (Figure  4),  for  different  sphere 
materials  (Figure  5), -for  different  azimuths,  6  (Figure  6),  and  for  different  sphere 
radii  (Figure  7).  The  degree  of  media  contrast  and  the  asymmetry  of  the  source 
determined  the  size  of  S  wave  generation.  For  the  most  extreme  case  investigated, 
we  found  that  the  maximum  S  wave  amplitude  in  the  time  domain  is  almost  3 
times  the  maximum  P  wave  amplitude  and  is  more  than  50%  greater  thaD  the  P 


wave  amplitude  for  the  center  explosion.  Before  extrapolating  these  time-domain 
results  to  more  realistic  environments,  one  must  remember  that  this  is  a  very  high 
frequency  source  with  a  contribution  to  the  far-field  displacement,  which  is  essen¬ 
tially  fiat  in  frequency. 

In  the  frequency  domain  the  asymmetry  of  the  source  region  has  significant 
effects  on  the  radiation  patterns  of  short  period,  but  has  only  minor  effects  on 
that  of  long  period.  We  also  found  that  the  the  maximum  long  period  P  wave 
amplitude  was  at  least  a  factor  of  two  greater  than  the  5  wave  maximum  ampli¬ 
tude  for  our  models.  We  also  investigated  the  amount  of  moment  error  caused  by 
using  a  finite  time  window.  With  care  one  can  obtain  meaningful  estimates  using 
only  a  short  duration  which  contains  relatively  few  pulses.  Since  the  moment 
relation,  which  is  obtained  in  the  limit  of  zero  frequency,  is  independent  of  the 
source  location  in  the  sphere  and  does  not  involve  the  geometry  or  dimensions  of 
the  sphere,  we  infer  that  it  probably  holds  for  embedded  homogeneous  source 
regions  of  any  shape. 
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APPENDIX  A 


In  this  appendix,  we  give  the  derivation  of  equations  (10)  and  (11) 
from  equations  (7 ),  (8)  and  (9).  The  addition  theorem  for 
zonal  harmonics  is 


i  (l— -mV. 

Pl (cos<5)-«  PAcoe^)^,/m(co8^oXco®m ^cosm ^Q+sinm $sinm ^qXAI) 


Substituting  (Al)  into  equation  (8),  we  obtain  equation  (10). 
Defining 

^n“APn^COS^e'T 

(A 2) 

,we  have 

— Pn  r0Xe  r 

(A3) 

Thus 

<t>n  i— Pn  (cos5)(sin0osin^ocos0— sinfein  ^cos0o) 
<t>n2’m~f> n  (cos^cos^osinSoos^— sin&ocos^>DCos0) 

4>n  a"'—?*  (cos5)(sin%sinfein(^— <^o) 

where 

a  coso 

Obviously, 

dPn  (cos  S)  , 

d<j>Q  ~P*  (cos^io^in0osin(^-0o)— 3, 

(A4) 

,  .  ,  (cos5) 

0Bl-sin0o  cote0cos<f>Q<f>n3} 

(A5) 

and 

dPn  (cos  6) 

<P  f»2"“cos^o  qq  ^ote^\n(t>Q<t>%3. 

(A6) 

From  equations  (Al)  -  (A6),  it  is  not  difficult  to  get  equation 

(11). 


APPENDIX  B 


In  this  appendix,  we  give  the  analytic  solution  for  the  source  at 
the  center  of  the  cavity. 

From  equation  (1),  we  have 

u  1-AL070(tair)-6L0+o(*air) 
u2“<*lo7o(*V) 

where  A  is  ik^g(u>)A  0.  Other  terms  in  (1)  are  zero.  From 
the  boundary  condition  at  r— a,  we  have 


I&1  a  ^o,3  (£2) 
l‘,JA  -yF  0+3(?i 


-*42),(«2)  *i2)'(fi) 

)  Jo(f.)  TFfcK,) 


where  £,  ■= A:Q  a ,  7* 


Vi**, 

Ma,’ 


According  to  the  relations  of  Bessel  functions,  we  have 

From  equations  (B2)  and  (B3),  u  2  can  be  written  as 


U2-w(w)S7<9r)2‘F<*“-r)e’' 
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Table  1.  Model 


a 

/? 

P 

radius 

(km  /sec ) 

(km  /sec) 

(g/cm3) 

km 

1.80 

0.41 

1.84 

0.5 

4.55 

2.57 

2.45 

oo 

* 
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FIGURE  CAPTIONS 


Figure  1.  The  geometry  and  the  coordinate  system  used  in  this  study,  star  is 
source,  and  triangle  is  receiver. 

Figure  2.  Comparisons  of  numerical  results  with  two  known  analytical  solutions. 
The  upper  is  for  the  wavefields  of  whole  space,  while  the  lower  is  for  the  wav. 
fields  from  a  centered  source  in  an  embedded  sphere.  '*  Num.  "  denotes  nu¬ 
merical  and  "  Ana.  "  denotes  analytic. 

Figure  3.  A  example  of  wavefields  from  an  off-centered  source  in  a  embedded 
sphere  with  a  radium  of  0.5  km.  "  P  "  is  for  the  P  waves,  "  S  "  is  for  the  S 
waves.  "  Rad  "  is  for  the  er  component,  "  Az  "  is  for  the  component.  The 
numbers  are  the  peak  amplitudes. 

Figure  4.  Wave  fields  for  different  source  distances  to  the  center  of  the  embedded 
sphere  with  a  radium  of  0.5  km.  The  10-2  on  the  second  and  third  traces  of 
the  first  line,  apply  to  all  the  traces  in  the  second  and  third  columns.  r0  is  the 
distance  from  the  source  to  the  center.  The  vertical  scale  is  same  for  all 
seismograms. 

Figure  5.  Wave  fields  for  different  materials  in  the  embedded  sphere,  a  is 
compressional  velocity,  and  /?  is  shear  velocity.  Nafe’s  velocity-density  rela¬ 
tion  is  assumed  for  the  shear  velocity  and  density  p. 

Figure  6.  Wave  fields  as  a  function  of  6  for  the  standard  model.  6  ranges  from  0 
to  90°. 

Figure  7.  Wave  fields-for  different  radii,  0.2  -  0.9  km,  of  the  embedded  sphere. 

Figure  8.  Radial  components  of  P  waves  at  and  180°  for  the  different  source 
locations.  In  the  cartoon  of  the  sphere  at  the  top,  the  star  indicates  the 
source. 

Figure  9.  Radial  components  of  P  waves  at  &-0  and  180°  for  different  P  and  5 
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velocities  and  density  in  the  sphere. 


Figure  10.  Spectral  radiation  patterns  at  periods  of  1  and  20  seconds  for  different 
source  locations.  "P"  is  for  the  radial  component  of  P  waves,  "S"  is  for  the  6 
component  of  S  waves.  On  the  cartoon  of  the  sphere  at  the  top,  the  direc¬ 
tion  of  6  is  indicated.  The  star  indicates  the  source.  "+"  is  the  center.  Dash 
circles  with  the  same  centers  as  the  P  radiation  patterns  are  the  P  wave  ra¬ 
diation  patterns  for  a  source  located  at  the  center.  The  numbers  on  the  right 
of  the  S  wave  radiation  patterns  are  the  multicative  factors  used  to  adjust 
the  the  size  of  S  wave  radiation  patterns  for  comparison  with  the  correspond¬ 
ing  P  wave  radiation  patterns  in  the  figure. 

Figure  11.  Spectral  radiation  patterns  for  different  elastic  media  of  the  embedded 
sphere. 

Figure  12.  Specral  radiation  patterns  for  different  radii  of  the  elastic  sphere. 

Figure  13.  Moment  estimates  as  a  function  of  record  length.  The  numbers  are  the 
percentage  of  estimated  moment  to  input  moment  obtained  using  the  time 
length  from  the  beginning  to  the  dash  line  mark  on  the  record  below  the 
number.  The  upper  trace  is  the  seismogram  of  a  centered  explosion.  The 
lower  one  is  the  seismogram  of  an  off-center  explosion.  Detailed  parameters 
are  given  in  the  text. 

Figure  14.  Ratio  of  the  observed  moment  to  the  input  source  moment,  Moit/M0, 
as  a  function  of  the  compressional  velocities  Oj  of  the  embedded  sphere  for 
different  whole-space  velocities,  Oj,  assuming  Nafe’s  velocity-density  relation. 
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Lexington,  MA  02173-0073  (3  copies) 
Prof  Fred  K.  Lamb 

University  of  Illinois  at  Urbana-Champaign 
Department  of  Physics 
1 1 10  West  Green  Street 
Urbana,  EL  61801 

Prof.  Charles  A.  Langston 
Geosciences  Department 
403  Deike  Building 
The  Pennsylvania  State  University 
University  Park,  PA  16802 

Prof.  Thome  Lay 

Institute  of  Tectonics 

Earth  Science  Board 

University  of  California,  Santa  Cruz 

Santa  Cruz,  CA  95064 

Prof.  Arthur  Lemer-Lam 
Lamont-Doherty  Geological  Observatory 
of  Columbia  University 
Palisades,  NY  10964 


Dr.  Christopher  Lynnes 
Teledyne  Geotech 
314  Montgomery  Street 
Alexandria,  V A  22314 


Prof.  Peter  Malin 
Department  of  Geology 
Old  Chemistry  Bldg. 

Duke  University 
Durham,  NC  27706 

Dr.  Randolph  Martin,  III 
New  England  Research,  Inc. 

76  Olcott  Drive 

White  River  Junction,  VT  05001 


Prof.  Thomas  V.  McEvilly 
Seismographic  Station 
University  of  California 
Berkeley,  CA  94720 


Dr.  Keith  L.  McLaughlin 
S-CUBED 

A  Division  of  Maxwell  Laboratory 
P.O.  Box  1620 
La  Jolla,  CA  92038-1620 
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Prof.  William  Menke 
Lamont-Doherty  Geological  Observatory 
of  Columbia  University 
Palisades,  NY  10964 


Stephen  Miller 
SRI  International 
333  Ravenswood  Avenue 
Box  AF  1 16 

Menlo  Park,  CA  94025-3493 

Prof.  Bernard  Minster 
IGPP,  A-025 

Scripps  Institute  of  Oceanography 
University  of  California,  San  Diego 
La  Jolla,  CA  92093 

Prof.  Brian  J.  Mitchell 

Department  of  Earth  &  Atmospheric  Sciences 
St.  Louis  University 
St.  Louis,  MO  63156 


Prof.  Charles  G.  Sammis 
Center  for  Earth  Sciences 
University  of  Southern  California 
University  Park 
Los  Angeles,  CA  90089-0741 

Prof.  Christopher  H.  Scholz 
Lamont-Doherty  Geological  Observatory 
of  Columbia  University 
Palisades,  NY  10964 


Thomas  J.  Sereno,  Jr. 

Science  Application  Int'l  Corp. 
10260  Campus  Point  Drive 
San  Diego,  CA  92121 


Prof.  David  G.  Simpson 
Lamont-Doherty  Geological  Observatory 
of  Columbia  University 
Palisades,  NY  10964 


Mr.  Jack  Murphy 

S-CUBED,  A  Division  of  Maxwell  Laboratory 
11800  Sunrise  Valley  Drive 
Suite  1212 

Reston,  VA  22091  (2  copies) 

Prof.  John  A.  Orcutt 
IGPP,  A-025 

Scripps  Institute  of  Oceanography 
University  of  California,  San  Diego 
La  Jolla,  CA  92093 

Prof.  Keith  Priestley 
University  of  Cambridge 
Bullard  Labs,  Dept,  of  Earth  Sciences 
Madingley  Rise,  Madingley  Rd. 

Cambridge  CB3  OEZ,  ENGLAND 

Dr.  Jay  J.  Pulli 
Radix  Systems,  Inc. 

2  Taft  Court,  Suite  203 
Rockville,  MD  20850 


Prof.  Paul  G.  Richards 
Lam  on  t  Doherty  Geological  Observatory 
of  Columbia  University 
Palisades,  NY  10964 


Dr.  Wilmer  Rivers 
Teledyne  Geotech 
314  Montgomery  Street 
Alexandria,  VA  22314 


Dr.  Jeffrey  Stevens 
S-CUBED 

A  Division  of  Maxwell  Laboratory 
P.O.  Box  1620 
La  Jolla,  CA  92038-1620 

Prof.  Brian  Stump 

Institute  for  the  Study  of  Earth  &  Man 
Geophysical  Laboratory 
Southern  Methodist  University 
Dallas,  TX  75275 

Prof.  Jeremiah  Sullivan 

University  of  Illinois  at  Urbana-Champaign 

Department  of  Physics 

1110  West  Green  Street 

Urbana,  EL  61801 

Prof.  Clifford  Thurber 
University  of  Wisconsin-Madison 
Department  of  Geology  &  Geophysics 
1215  West  Dayton  Street 
Madison,  WS  53706 

Prof.  M.  Nafi  Toksoz 
Earth  Resources  Lab 
Massachusetts  Institute  of  Technology 
42  Carleton  Street 
Cambridge,  MA  02142 

Prof.  John  E.  Vidale 

University  of  California  at  Santa  Cruz 

Seismological  Laboratory 

Santa  Cruz,  CA  95064 
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Prof.  Terry  C.  Wallace 
Department  of  Geosciences 
Building  #77 
University  of  Arizona 
Tucson,  AZ  85721 

Dr.  William  Wortman 
Mission  Research  Corporation 
8560  Cinderbed  Rd. 

Suite  #  700 

Newington,  VA  22122 

Prof.  Francis  T.  Wu 
Department  of  Geological  Sciences 
State  University  of  New  York 
it  Binghamton 
Vestaf,  NY  1 390 1 


OTHERS  (United  States) 


Dr.  Monem  Abdel-Gawad 
Rockwell  International  Science  Center 
1049  Camino  Dos  Rios 
Thousand  Oaks,  CA  91360 


Prof.  Keiiti  Aid 

Center  for  Earth  Sciences 

University  of  Southern  California 

University  Park 

Los  Angeles,  CA  90089-0741 

Prof.  Shelton  S.  Alexander 
Geosciences  Department 
403  Deike  Building 
The  Pennsylvania  State  University 
University  Park,  PA  16802 

Dr.  Kenneth  Anderson 
BBNSTC 
Mail  Stop  14/IB 
Cambridge,  MA  02238 


Dr.  Ralph  Archuleta 
Department  of  Geological  Sciences 
University  of  California  at  Santa  Barbara 
Santa  Barbara,  CA  93102 


Dr.  Susan  Beck 
Department  of  Geosciences 
Bldg.  #  77 

University  of  Arizona 
Tucson,  AZ  85721 

Dr.  T.J.  Bennett 
S-CUBED 

A  Division  of  Maxwell  Laboratory 
1 1800  Sunrise  Valley  Drive,  Suite  1212 
Reston,  VA  22091 

Mr.  William  J.  Best 
907  Westwood  Drive 
Vienna,  VA  22180 


Dr.  N.  Biswas 
Geophysical  Institute 
University  of  Alaska 
Fairbanks,  AK  99701 


Dr.  Stephen  Bratt 
Center  for  Seismic  Studies 
1300  North  17th  Street 
Suite  1450 

Arlington,  VA  22209 


Michael  Browne 
Teledyne  Geotech 
3401  Shiloh  Road 
Garland,  TX  75041 


Mr.  Roy  Burger 
1221  Serry  Road 
Schenectady,  NY  12309 


Dr.  Robert  Burridge 
Schlumberger-Doll  Research  Center 
Old  Quarry  Road 
Ridgefield,  CT  06877 


Dr.  W.  Winston  Chan 
Teledyne  Geotech 
314  Montgomery  Street 
Alexandria,  VA  22314-1581 


Dr.  Theodore  Cherry 
Science  Horizons,  Inc. 

710  Encinitas  Blvd.,  Suite  200 
Encinitas,  CA  92024  (2  copies) 


Prof.  Jon  F.  Claerbout 
Department  of  Geophysics 
Stanford  University 
Stanford,  CA  94305 


Prof.  Robert  W.  Clayton 

Seismological  Laboratory 

Division  of  Geological  &  Planetary  Sciences 

California  Institute  of  Technology 

Pasadena,  CA  91125 

Prof.  F.  A.  Dahlen 

Geological  and  Geophysical  Sciences 

Princeton  University 

Princeton,  NJ  08544-0636 


Mr.  Charles  Doll 
Earth  Resources  Laboratory 
Massachusetts  Institute  of  Technology 
42  Carleton  St. 

Cambridge,  MA  02142 

Prof.  Adam  Dziewonski 
Hoffman  Laboratory,  Harvard  Univ. 

Dept,  of  Earth  Atmos.  &  Planetary  Sciences 
20  Oxford  St 
Cambridge,  MA  02138 
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Prof.  John  Ebel 

Department  of  Geology  &  Geophysics 
Boston  College 
Chestnut  Hill,  MA  02167 


Eric  Fielding 
SNEE  Hall 
1NSTOC 

Cornell  University 
Ithaca,  NY  14853 

Dr.  John  Foley 

Phillips  Laboratory/LWH 

Hanscom  AFB,  MA  01731-5000 


Prof.  Donald  Forsyth 
Department  of  Geological  Sciences 
Brown  University 
Providence,  RI  02912 


Dr.  Cliff  Frolich 
Institute  of  Geophysics 
8701  North  Mopac 
Austin,  TX  78759 


Dr.  Anthony  Gangi 
Texas  A&M  University 
Department  of  Geophysics 
College  Station,  TX  77843 


Dr.  Freeman  Gilbert 
IGPP,  A-025 

Scripps  Institute  of  Oceanography 
University  of  California 
La  Jolla,  CA  92093 

Mr.  Edward  Giller 
Pacific  Sierra  Research  Corp. 

1401  Wilson  Boulevard 
Arlington,  VA  22209 


Dr.  Jeffrey  W.  Given 
SAIC 

10260  Campus  Point  Drive 
San  Diego,  CA  92121 


Prof.  Stephen  Grand 
University  of  Texas  at  Austin 
Department  of  Geological  Sciences 
Austin,  TX  78713-7909 


Prof.  Roy  Greenfield 
Geosciences  Department 
403  Deike  Building 
The  Pennsylvania  State  University 
University  Park,  PA  16802 

Dan  N.  Hagedorn 
Battelle 

Pacific  Northwest  Laboratories 
Battelle  Boulevard 
Richland,  WA  99352 

Dr.  James  Hannon 

Lawrence  Livermore  National  Laboratory 
P.  O.  Box  808 
Livermore,  CA  94550 


Prof.  Robert  B.  Herrmann 

Dept,  of  Earth  &  Atmospheric  Sciences 

St.  Louis  University 

St.  Louis,  MO  63156 


Ms.  Heidi  Houston 
Seismological  Laboratory 
University  of  California 
Santa  Cruz,  CA  95064 


Kevin  Hutchenson 
Department  of  Earth  Sciences 
St.  Louis  University 
3507  Laclede 
St.  Louis,  MO  63103 

Dr.  Hans  Israelsson 
Center  for  Seismic  Studies 
1300  N.  17th  Street,  Suite  1450 
Arlington,  VA  22209-2308 


Prof.  Thomas  H.  Jordan 
Department  of  Earth,  Atmospheric 
and  Planetary  Sciences 
Massachusetts  Institute  of  Technology 
Cambridge,  MA  02139 

Prof.  Alan  Kafka 

Department  of  Geology  &  Geophysics 
Boston  College 
Chestnut  Hill,  MA  02167 


Robert  C.  Kemerait 
ENSCO,  Inc. 

445  Pineda  Court 
Melbourne,  FL  32940 


William  Kikendall 
Teledyne  Geotech 
3401  Shiloh  Road 
Garland,  TX  75041 


Prof.  Amos  Nur 
Department  of  Geophysics 
Stanford  University 
Stanford,  CA  94305 


Prof.  Leon  Knopoff 

University  of  California 

Institute  of  Geophysics  &  Planetary  Physics 

Los  Angeles,  CA  90024 


Prof.  Jack  Oliver 
Department  of  Geology 
Cornell  University 
Ithaca,  NY  14850 


Prof.  John  Kuo 

Aldridge  Laboratory  of  Applied  Geophysics 
Columbia  University 
842  Mudd  Bldg. 

New  York,  NY  10027 

Prof.  L.  Tirr-'thy  Long 
School  of  Geophysical  Sciences 
Georgia  Institute  of  Technology 
Atlanta,  GA  30332 


Dr.  Gary  McCartor 
Department  of  Physics 
Southern  Methodist  University 
Dallas,  TX  75275 


Prof.  An  McGarr 
Mail  Stop  977 
Geological  Survey 
345  Middlefield  Rd. 
Menlo  Park,  CA  94025 

Dr.  George  Mellman 
Sierra  GeoDhysics 
11255  Kirkland  Way 
Kirkland,  WA  98033 


Prof.  John  Nabelek 
College  of  Oceanography 
Oregon  State  University 
Corvallis,  OR  97331 


Prof.  Geza  Nagy 

University  of  California,  San  Diego 
Department  of  Ames,  M  S.  B-010 
La  Jolla,  CA  92093 


Dr  Keith  K.  Nakanishi 

Lawrence  Livermore  National  Laboratory 

L-205 

P.  O.  Box  808 
Livermore,  CA  94550 


Dr.  Kenneth  Olsen 
P.  O.  Box  1273 
Linwood,  WA  98046-1273 


Prof.  Jeffrey  Park 

Department  of  Geology  and  Geophysics 
Kline  Geology  Laboratory 
P.  O.  Box  6666 
New  Haven,  CT  065 1 1  -8 1 30 

Howard  J.  Patton 

Lawrence  Livermore  National  Laboratory 
L-205 

P.  O.  Box  808 
Livermore,  CA  94550 

Prof.  Robert  Phinney 
Geological  &  Geophysical  Sciences 
Princeton  University 
Princeton,  NJ  08544-0636 


Dr.  Paul  Pomeroy 
Rondout  Associates 
P.O.  Box  224 
Stone  Ridge,  NY  12484 


Dr.  Norton  Rimer 
S-CUBED 

A  Division  of  Maxwell  Laboratory 
P.O.  Box  1620 
La  Jolla,  C A  92038-1620 

Prof.  Larry  J.  Ruff 
Department  of  Geological  Sciences 
1006  C.C.  Little  Building 
University  of  Michigan 
Ann  Arbor,  MI  48109-1063 

Dr.  Richard  Sailor 
TASC  Inc. 

55  Walkers  Brook  Drive 
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Dr.  Susan  Schwartz 

Institute  of  Tectonics 

1156  High  St. 

Santa  Cruz,  CA  95064 

Dr.  David  Taylor 

ENSCO,  Inc. 

445  Pineda  Court 

Melbourne,  FL  32940 

John  Sherwin 

Teledyne  Geotech 

3401  Shiloh  Road 

Garland,  TX  75041 

Dr.  Steven  R.  Taylor 

Lawrence  Livermore  National  Laboratory 

L-205 

P.  O.  Box  808 

Livermore,  CA  94550 

Dr.  Matthew  Sibol 

Virginia  Tech 

Seismological  Observatory 

4044  Derring  Hall 

Blacksburg,  VA  24061-0420 

Professor  Ta-Liang  Teng 

Center  for  Earth  Sciences 

University  of  Southern  California 

University  Park 

Los  Angeles,  CA  90089-0741 

Dr.  Albert  Smith 

Lawrence  Livermore  National  Laboratory 

L-205 

P.  O.  Box  808 

Livermore,  CA  94550 

Dr.  Gregory  van  der  Vink 

IRIS,  Inc. 

1616  North  Fort  Myer  Drive 

Suite  1440 

Arlington,  VA  22209 

Prof.  Robert  Smith 

Department  of  Geophysics 

University  of  Utah 

1400  East  2nd  South 

Salt  Lake  City,  UT  84112 

Professor  Daniel  Walker 

University  of  Hawaii 

Institute  of  Geophysics 

Honolulu,  HI  96822 

Dr.  Stewart  W.  Smith 

Geophysics  AK-50 

University  of  Washington 

Seattle,  WA  98195 

William  R.  Walter 

Seismological  Laboratory 

University  of  Nevada 

Reno,  NV  89557 

Donald  L.  Springer 

Lawrence  Livermore  National  Laboratory 

L-205 

P.  O.  Box  808 

Livermore,  CA  94550 

Dr.  Raymond  Willeman 

Phillips  Laboratory  -  OL-AA/LWH 

Hanscom  AFB,  MA  01731-5000 

Dr.  George  Sutton 

Rondout  Associates 

P.O.  Box  224 

Stone  Ridge,  NY  12484 

Dr.  Gregory  Wojcik 

Weidlinger  Associates 

4410  El  Camino  Real 

Suite  1 10 

Los  Altos,  CA  94022 

Prof.  L.  Sykes 

Lamont-Doherty  Geological  Observatory 
of  Columbia  University 

Palisades,  NY  10964 

Dr.  Lorraine  Wolf 

Phillips  Laboratory/LWH 

Hanscom  AFB,  MA  01731-5000 

Prof.  Pradeep  Talwani 

Department  of  Geological  Sciences 

University  of  South  Carolina 

Columbia,  SC  29208 

Dr.  Gregory  B.  Young 

ENSCO,  Inc. 

5400  Port  Royal  Road 

Springfield,  VA  22151-2388 
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Dr.  Eileen  Vergino 

Lawrence  Livermore  National  Laboratory 
L-205 

P.  O.  Box  808 
Livermore,  CA  94550 

J.  J.  Zucca 

Lawrence  Livermore  National  Laboratory 
P.  O.  Box  808 
Livermore,  CA  94550 


GOVERNMENT 


Dr.  Ralph  Alewine  III 

DARPA/NMRO 

1400  Wilson  Boulevard 

Arlington,  VA  22209-2308 

Dr.  Dale  Glover 

DIA/DT-1B 

Washington,  DC  20301 

Mr.  Janies  C.  Battis 

Phillips  Laboratory/LWH 

Hanscom  AFB,  MA  01731-5000 

Dr.  T.  Hanks 

USGS 

Nat'l  Earthquake  Research  Center 
345  Middlefield  Road 

Menlo  Park,  CA  94025 

Harley  Benz 

U.S.  Geological  Survey,  MS-977 

345  Middlefield  Rd. 

Menlo  Park,  CA  94025 

Dr.  Roger  Hansen 

AFTAGTT 

Patrick  AFB,  FL  32925 

Dr.  Robert  Blandford 

AFTACAT 

Center  for  Seismic  Studies 

1300  North  17th  St.  Suite  1450 

Arlington,  VA  22209-2308 

Paul  Johnson 

ESS-4,  Mail  Stop  J979 

Los  Alamos  National  Laboratory 

Los  Alamos,  NM  87545 

Eric  Chael 

Division  9241 

Sandia  Laboratory 

Albuquerque,  NM  87185 

Janet  Johnston 

Phillips  Laboratory/LWH 

Hanscom  AFB,  MA  01731-5000 

Dr.  John  J.  Cipar 

Phillips  Laboratory/LWH 

Hanscom  AFB,  MA  01731-5000 

Dr.  Katharine  Kadinsky-Cade 
Phillips  Laboratory/LWH 

Hanscom  AFB,  MA  01731-5000 

Cecil  Davis 

Group  P-15,  Mail  Stop  D406 

P.O.  Box  1663 

Los  Alamos  National  Laboratory 

Los  Alamos,  NM  87544 

Ms.  Ann  Kerr 

IGPP,  A-025 

Scripps  Institute  of  Oceanography 
University  of  California,  San  Diego 
La  Jolla,  CA  92093 

Mr.  Jeff  Duncan 

Office  of  Congressman  Markey 

2133  Rayburn  House  Bldg. 

Washington,  DC  20515 

Dr.  Max  Koontz 

US  Dept  of  Energy/DP  5 

Forrestal  Building 

1000  Independence  Avenue 
Washington,  DC  20585 

Dr.  Jack  Evemden 

USGS  -  Earthquake  Studies 

345  Middlefield  Road 

Menlo  Park,  CA  94025 

Dr.  W.H.K.  Lee 

Office  of  Earthquakes,  Volcanoes, 
&  Engineering 

345  Middlefield  Road 

Menlo  Park,  CA  94025 

Art  Frankel 

USGS 

922  National  Center 

Reston,  VA  22092 

Dr.  William  Leith 

U.S.  Geological  Survey 

Mail  Stop  Q28 

Reston,  VA  22092 
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Dr.  Richard  Lewis  Mr.  Chris  Paine 

Director,  Earthquake  Engineering  &  Geophysics  Office  of  Senator  Kennedy 

U.S.  Army  Corps  of  Engineers  SR  315 

Box  63 1  United  States  Senate 

Vicksburg,  MS  39180  Washington,  DC  20510 


James  F.  Lewkowicz 
Phillips  Laboratory/LWH 
Hanscom  AFB,  MA  01731-5000 


Colonel  Jerry  J.  Perrizo 
AFOSR/NP,  Building  410 
Bolling  AFB 

Washington,  DC  20332-6448 


Mr.  Alfred  Lieberman 
ACDA/VI-OA'State  Department  Bldg 
Room  5726 
320-  21st  Street,  NW 
Washington,  DC  20451 

Stephen  Mangino 
Phillips  Laboratory/LWH 
Hanscom  AFB,  MA  01731-5000 


Dr.  Robert  Masse 
Box  25046,  Mail  Stop  967 
Denver  Federal  Center 
Denver,  CO  80225 


Art  McGarr 

U.S.  Geological  Survey,  MS-977 
345  Middlefield  Road 
Menlo  Park,  CA  94025 


Richard  Morrow 
ACDA/VI,  Room  5741 
320  21st  Street  N.W 
Washington,  DC  20451 


Dr.  Carl  Newton 

Los  Alamos  National  Laboratory 

P.O.  Box  1663 

Mail  Stop  C335,  Group  ESS-3 
Los  Alamos,  NM  87545 

Dr.  Bao  Nguyen 

AFTAC/TTR 

Patrick  AFB,  FL  32925 


Dr.  Kenneth  H.  Olsen 

Los  Alamos  Scientific  Laboratory 

P.  O.  Box  1663 

Mail  Stop  D-406 

Los  Alamos,  NM  87545 


Dr.  Frank  F.  Pilotte 

HQ  AFTAC/TT 

Patrick  AFB,  FL  32925-6001 


Katie  Poley 
CLA-ACIS/TMC 
Room  4X16NHB 
Washington,  DC  20505 


Mr.  Jack  Rachlin 
U.S.  Geological  Survey 
Geology,  Rm  3  Cl 36 
Mail  Stop  928  National  Center 
Reston,  VA  22092 

Dr.  Robert  Reinke 
WL/NTESG 

Kirtland  AFB,  NM  87117-6008 


Dr.  Byron  Ristvet 

HQ  DNA,  Nevada  Operations  Office 

Attn:  NVCG 

P.O.  Box  98539 

Las  Vegas,  NV  89193 

Dr.  George  Rothe 
HQ  AFTAC/TTR 
Patrick  AFB,  FL  32925-6001 


Dr.  Alan  S.  Rvall,  Jr. 
DARPA/NMRO 
1400  Wilson  Boulevard 
Arlington,  VA  22209-2308 


Dr.  Michael  Shore 
Defense  Nuclear  Agency/SPSS 
6801  Telegraph  Road 
Alexandria,  VA  22310 


it 


Mr.  Charles  L.  Taylor 
Phillips  Laboratory/LWH 
Hanscom  AFB,  MA  01731-5000 


Dr.  Larry  Turnbull 
CIA-OSWR/NED 
Washington,  DC  20505 


Dr.  Thomas  Weaver 
Los  Alamos  National  Laboratory 
P.O.  Box  1663,  Mail  Stop  C335 
Los  Alamos,  NM  87545 


Phillips  Laboratory 
Research  Library 
ATTN:  SULL 

Hanscom  AFB  ,  MA  01731-5000  (2  copies) 


Phillips  Laboratory 
Attn:  XO 

Hanscom  AFB,  MA  01731-5000 


Phillips  Laboratory 
Attn:  LW 

Hanscom  AFB,  MA  01731-5000 


DARPA/PM 

1400  Wilson  Boulevard 

Arlington,  VA  22209 


Defense  Technical  Information  Center 
Cameron  Station 

Alexandria,  V A  223 14  (5  copies) 


Phillips  Laboratory 
ATTN:  SUL 

Kirtland  AFB,  NM  87 1 17-6008 


Defense  Intelligence  Agency 

Directorate  for  Scientific  &  Technical  Intelligence 

Atm:  DT1B 

Washington,  DC  20340-6158 


Secretary  of  the  Air  Force 
(SAFRD) 

Washington,  DC  20330 


AFTAC/CA 

(STINFO) 

Patrick  AFB,  FL  32925-6001 


TACTEC 

Batteile  Memorial  Institute 
505  King  Avenue 

Columbus,  OH  43201  (Final  Report  Only) 


HQ  DNA 

Attn:  Technical  Library 
Washington,  DC  20305 


Office  of  the  Secretary  Defense 
DDR  &  E 

Washington,  DC  20330 


D  AR  P  A/RMO/RETRIE  V  AL 
1400  Wilson  Boulevard 
Arlington,  VA  22209 


DARPA/RMO/Security  Office 
1400  Wilson  Boulevard 
Arlington,  VA  22209 
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CONTRACTORS  (FOREIGN) 


Dr.  Ramon  Cabre,  S  J. 
Observatorio  San  Calixto 
Casilla  5939 
La  Paz,  Bolivia 


Prof.  Hans-Peter  Harjes 
Institute  for  Geophysik 
Ruhr  University /Bochum 
P.O.  Box  102148 
4630  Bochum  1,  FRG 

Prof.  Eystein  Husebye 
NTNF/NORSAR 
P.O.  Box  51 

N-2007  Kjeller,  NORWAY 


Prof.  Brian  L.N.  Kennett 
Research  School  of  Earth  Sciences 
Institute  of  Advanced  Studies 
G.P.O.  Box  4 

Canberra  2601,  AUSTRALIA 

Dr.  Bernard  Massinon 

Societe  Radiomana 

27  rue  Claude  Bernard 

75005  Paris,  FRANCE  (2  Copies) 


Dr.  Pierre  Mecheler 
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